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Abstract 

In  this  thesis  we  present  LaSLAT,  a  sensor  network  algorithm  that  uses 
range  measurements  between  sensors  and  a  moving  target  to  simulta¬ 
neously  localize  the  sensors,  calibrate  sensing  hardware,  and  recover 
the  target’s  trajectory. 

LaSLAT  is  based  on  a  Bayesian  filter  that  updates  a  probability 
distribution  over  the  parameters  of  interest  as  measurements  arrive. 
The  algorithm  is  distributable  and  requires  a  fixed  amount  of  storage 
space  with  respect  to  the  number  of  measurements  it  has  incorporated. 
LaSLAT  is  easy  to  adapt  to  new  types  of  hardware  and  new  physical 
environments  due  to  its  use  of  intuitive  probability  distributions:  one 
adaptation  demonstrated  in  this  thesis  uses  a  mixture  measurement 
model  to  detect  and  compensate  for  bad  acoustic  range  measurements 
due  to  echoes. 

We  present  results  from  a  centralized  implementation  of  LaSLAT 
using  a  network  of  Cricket  sensors.  In  both  2D  and  3D  networks, 
LaSLAT  is  able  to  localize  sensors  to  within  several  centimeters  of  their 
ground  truth  positions  while  recovering  a  range  measurement  bias  for 
each  sensor  and  the  complete  trajectory  of  the  mobile. 

Thesis  Supervisor:  Jonathan  Bachracli 
Title:  Research  Scientist 
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Chapter  1 

Introduction 


This  thesis  presents  LaSLAT,  an  algorithm  for  localizing  and  calibrat¬ 
ing  a  sensor  network  while  simultaneously  tracking  a  moving  target. 
LaSLAT  performs  these  tasks  using  range  measurements  between  the 
sensors  and  the  target.  LaSLAT  does  not  require  sensor  nodes  with 
known  positions,  ranging  information  between  sensors,  or  constraints 
on  the  path  of  the  target. 

1.1  Localization 

Many  sensor  network  applications  require  that  the  locations  of  the  in¬ 
dividual  sensors  be  known,  since  sensor  readings  are  in  general  of  little 
use  without  geographic  context.  However,  the  same  attributes  that 
make  sensor  networks  attractive  make  obtaining  this  information  diffi¬ 
cult.  By  placing  a  large  number  of  relatively  cheap  sensors,  it  is  possible 
to  obtain  many  accurate  measurements  from  sensors  close  to  phenom¬ 
ena  of  interest;  however,  the  sheer  number  of  sensors  and  the  need  to 
minimize  costs  precludes  manually  recording  the  sensors’  locations.  It 
also  precludes  brute  force  solutions  such  as  equipping  each  sensor  with 
a  GPS  unit. 

Consequently,  we  would  like  the  sensors  to  determine  their  own  po¬ 
sitions  after  placement.  This  is  known  as  localization,  and  is  typically 
achieved  by  having  each  sensor  compute  range  measurements  to  its 
neighboring  sensors,  then  algorithmically  embedding  the  graph  formed 
by  these  ranges  into  a  coordinate  system  (see  chapter  2).  This  coordi¬ 
nate  system  is  then  used  to  perform  location-dependent  tasks  such  as 
geographic  packet  routing  or  target  tracking. 
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Localization  is  often  complicated  by  the  difficulty  of  obtaining  enough 
accurate  pair-wise  range  measurements  between  sensors.  Inter-sensor 
ranges  can  be  corrupted  by  noise  or  lost  entirely  due  to  occluded  line- 
of-siglit.  Thus,  consistently  accurate  localization  requires  robustness  in 
the  face  of  missing  or  low  quality  measurements. 

Nevertheless,  localization  is  rarely  if  ever  the  purpose  of  a  network. 
Sensor  networks  are  typically  deployed  to  observe  active  phenomena  in 
the  environment,  and  require  accurate  localization  as  a  means  to  that 
end.  As  a  result,  there  is  pressure  in  localization  research  to  achieve 
accuracy  and  robustness  using  as  little  hardware  as  possible. 

1.2  Tracking 

Target  tracking  is  one  of  the  motivating  localization-dependent  appli¬ 
cations  of  sensor  networks.  In  tracking  applications,  sensors  jointly 
observe  phenomena,  which  may  be  people  or  objects  passing  through 
the  network  or  physical  effects  such  as  bullet  shock-waves  or  anoma¬ 
lous  sounds.  Once  a  phenomenon  is  detected,  the  sensors  collaborate 
to  determine  its  spatial  location.  This  estimate  is  reported  to  a  com¬ 
puter  or  person  monitoring  the  network.  Target  tracking  networks  can 
provide  indoor  navigation  services  to  hand-held  users  or  mobile  robots, 
track  friendlies  and  hostiles  on  a  battlefield,  or  monitor  movement  of 
inventory  in  a  store  or  warehouse. 

Tracking  is  a  well- understood  problem  (see  chapter  2).  Given  the  lo¬ 
cations  of  the  sensors  and  accurate  range  information  to  the  target,  it  is 
straightforward  to  determine  the  target’s  position.  Consequently,  tra¬ 
ditional  tracking  applications  tend  to  be  split  into  two  separate  phases. 
First,  the  network  is  localized  using  a  specialized  algorithm.  After 
localization  completes,  the  network  enters  a  target  tracking  phase  in 
which  target  positions  are  estimated  based  on  the  discovered  sensor 
positions. 


1.3  Calibration 

In  an  ideal  world,  sensors  would  arrive  from  the  factory  fully  calibrated 
to  begin  taking  accurate  measurements  of  their  surroundings.  However, 
this  ideal  situation  is  rarely  achieved.  For  instance,  deployment  con¬ 
ditions  such  as  temperature  affect  the  accuracy  of  ranging  algorithms 
based  on  acoustic  time-of-flight  by  altering  the  speed  of  sound.  Fur¬ 
thermore,  as  shown  in  [35],  differences  between  sensors  can  also  result 


in  mis-calibrations  that  are  difficult  to  correct  before  deployment.  Cali¬ 
bration  in  the  field  can  therefore  offer  meaningful  improvement  in  both 
localization  and  target  tracking  accuracy.  As  with  localization,  there  is 
considerable  economic  incentive  to  develop  auto-calibration  algorithms 
that  allow  sensors  to  self-calibrate  in  the  field  without  external  inter¬ 
vention. 


1.4  Simultaneous  Localization  and  Track¬ 
ing 

For  sensor  networks  deployed  to  track  moving  targets,  some  authors 
have  suggested  using  a  moving  target  (or  mobile)  to  assist  in  localizing 
the  network  [8,12,27].  This  approach  is  attractive  for  several  reasons: 

•  It  requires  no  additional  localization-specific  hardware  on  the  in¬ 
dividual  themselves,  potentially  reducing  both  their  size  and  cost. 

•  Unlike  sensors,  a  mobile  can  move  freely  in  the  network.  This  pro¬ 
vides  a  larger  number  and  a  greater  diversity  of  measurements  for 
use  in  localization,  which  helps  reduce  the  effect  of  noisy  measure¬ 
ments  and  static  environmental  obstacles. 

•  Use  of  a  mobile  means  that  the  sensors  are  no  longer  constrained 
to  have  line  of  sight  to  each  other.  This  enables  a  broader  range 
of  network  deployments:  in  particular,  it  allows  sensors  to  be 
placed  to  optimize  tracking  coverage  rather  than  to  facilitate  lo¬ 
calization.  This  is  especially  significant  because  many  types  of 
ranging  hardware  are  directional.  Traditional  localization  algo¬ 
rithms  therefore  introduce  a  deployment  conflict:  sensors  must 
be  oriented  so  that  they  can  range  to  several  other  sensors  in  ad¬ 
dition  to  observing  areas  where  targets  are  likely  to  appear.  This 
wastes  sensing  resources  and  limits  deployment  options.  SLAT 
eliminates  this  problem,  since  it  works  best  when  sensors  use 
their  full  field  of  view  to  sense  targets. 

•  Localization  using  a  mobile  potentially  allows  localization  esti¬ 
mates  to  continuously  improve,  even  after  the  network  begins 
tracking  targets.  In  particular,  sensors  that  are  closest  to  high 
traffic  areas  can  be  localized  very  precisely.  In  turn,  this  facili¬ 
tates  high  precision  target  tracking  in  these  areas. 

Existing  methods  fail  to  realize  the  full  potential  of  mobile-based 
localization.  In  particular,  most  require  that  the  position  of  the  mobile 
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be  known  at  all  times.  Consequently,  the  mobile  must  be  constrained  to 
a  well-known  trajectory  or  equipped  with  a  GPS-like  system.  Further¬ 
more,  existing  methods  ignore  the  sensor  calibration  problem  entirely. 

We  consider  a  more  general  scenario,  in  which  the  mobile  is  allowed 
to  move  arbitrarily  on  an  unknown  path.  As  it  moves,  the  mobile  pe¬ 
riodically  emits  a  signal  that  allows  sensors  in  the  network  to  compute 
their  range  to  the  mobile.  These  ranges  are  expected  to  be  corrupted 
by  noise.  Using  solely  these  noisy  range  measurements,  the  objective 
is  to  localize  the  sensors  in  the  network,  calibrate  the  sensors’  ranging 
hardware  if  necessary,  and  track  the  position  of  the  mobile.  Accord¬ 
ingly,  we  label  this  problem  Simultaneous  Localization  and  Tracking 
(SLAT). 

1.5  LaSLAT 

Our  solution  to  the  SLAT  problem  takes  the  form  of  a  Bayesian  filter. 
The  filter  uses  range  measurements  taken  by  the  network  to  update 
a  joint  probability  distribution  over  the  positions  of  the  sensors,  the 
trajectory  of  the  mobile,  and  the  calibration  parameters  of  the  network. 
To  avoid  some  of  the  representational  and  computational  complexity  of 
general  Bayesian  filtering,  we  use  Laplace’s  method  to  approximate  our 
state  with  a  Gaussian  after  incorporating  each  batch  of  measurements. 
As  a  result,  we  call  our  algorithm  LaSLAT. 

LaSLAT  inherits  many  desirable  properties  from  the  Bayesian  frame¬ 
work.  The  probabilistic  model  used  in  LaSLAT  insures  that  measure¬ 
ment  noise  is  averaged  out  as  more  measurements  become  available, 
which  improves  localization  and  tracking  accuracy  in  high-traffic  areas. 
The  filtering  framework  incorporates  measurements  in  small  batches, 
providing  on-line  estimates  of  all  locations,  calibration  parameters,  and 
their  uncertainties.  This  allows  the  mobile  to  be  tracked  in  near  real¬ 
time.  It  also  speeds  up  the  convergence  of  the  algorithm  and  reduces 
impact  on  the  network. 

Since  LaSLAT  maintains  estimate  uncertainties,  mobiles  can  be  dis¬ 
patched  on-line  if  desired  to  improve  localization  estimates  in  regions 
where  localization  uncertainty  is  high.  As  we  required  in  the  prob¬ 
lem  statement,  mobiles  may  move  arbitrarily  through  the  environment, 
with  no  constraint  on  their  trajectory  or  velocity.  Furthermore,  multi¬ 
ple  mobiles  may  be  used  in  conjunction  to  expedite  initial  localization. 

If  available,  ancillary  localization  information  such  as  position  es¬ 
timates  from  GPS,  anchor  nodes,  or  radio-based  ranging  can  be  eas¬ 
ily  incorporated  into  our  framework.  Our  algorithm  avoids  expensive 
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matrix  operations  by  operating  on  sparse  inverse  covariance  matrices 
rather  than  operating  directly  on  dense  covariance  matrices.  This  helps 
LaSLAT  run  quickly  and  facilitates  performing  the  LaSLAT  computa¬ 
tions  distributedly  in  the  network. 

The  user  deploying  the  network  may  specify  a  coordinate  system 
for  LaSLAT  to  honor  when  localizing  sensors.  If  no  coordinate  system 
is  specified,  LaSLAT  recovers  locations  in  a  coordinate  system  that  is 
correct  up  to  a  translation,  rotation,  and  possible  reflection. 

We  demonstrate  these  features  by  accurately  localizing  a  dense  net¬ 
work  of  27  sensors  to  within  one  or  two  centimeters.  The  sensor  nodes 
are  wireless  Crickets  [1]  capable  of  measuring  their  distance  to  a  mov¬ 
ing  beacon  using  a  combination  of  ultrasound  and  radio  pulses.  In  a 
larger  and  sparser  network,  we  localize  sensors  to  within  about  eight 
centimeters.  In  both  cases,  a  measurement  bias  parameter  is  accurately 
calibrated  for  all  nodes.  Finally,  we  present  results  from  two  experi¬ 
ments  in  three  dimensions,  in  which  nodes  were  localized  to  within 
seven  centimeters. 
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Chapter  2 

Related  Work 

2.1  Localization 

Localization  is  a  well  established  problem  in  sensor  networks,  and  many 
approaches  have  been  developed  in  the  literature.  The  vast  majority 
of  these  techniques  treat  localization  as  a  stand-alone  procedure  that 
takes  place  in  advance  of  or  concurrently  with  other  operations  in  the 
network.  In  the  traditional  localization  problem,  sensors  are  capable  of 
determining  the  distance  between  themselves  and  other  nearby  sensors. 
This  may  be  done  using  relatively  accurate  methods  such  as  acoustic 
time-of- flight  (TOF),  or  using  imprecise  but  cheap  methods  such  as 
radio  signal  strength  indicators  (RSSI).  In  dense  networks,  radio  hop 
count  (DV)  may  be  used  as  a  surrogate  for  distance. 

Several  authors  [10,17,31]  have  used  multidimensional  scaling  (MDS) 
as  a  localization  technique.  MDS  is  very  effective  at  recovering  sensor 
topologies  when  accurate  distances  are  available  between  all  pairs  of 
sensors.  Its  performance  degrades  substantially  when  some  distances 
are  unavailable. 

Many  algorithms  [4,5,20,23,25,30,32]  also  depend  on  the  presence  of 
anchor  nodes,  which  are  sensors  that  know  their  true  locations  a  priori 
through  out-of-band  means  such  as  pre-progrannning  or  GPS.  These 
algorithms  localize  non-anchor  nodes  by  construction  using  distances 
from  each  non-anchor  node  to  several  anchor  nodes.  However,  anchor 
nodes  are  expensive,  and  not  always  perfectly  accurate  (i.e.  GPS). 
Furthermore,  many  of  these  algorithms  suffer  from  propagation  of  er¬ 
rors:  sensors  several  hops  from  a  beacon  often  accumulate  considerable 
position  error. 

Priyantha  et  al.  [29]  suggest  a  two-phase  approach  that  is  similar 
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in  some  ways  to  the  approach  proposed  in  this  thesis.  An  initialization 
phase  uses  an  imprecise  technique  based  on  radio  connectivity  between 
nodes  to  embed  the  sensors  in  a  plane.  Once  initialization  is  complete,  a 
second  phase  uses  precise  time  of  flight  ranges  between  adjacent  sensors 
to  make  local  position  adjustments. 

Ihler  et  al.  [16]  treat  localization  as  an  inference  problem  on  a  graph¬ 
ical  model.  This  allows  them  to  use  nonparanretric  belief  propagation 
(NBP)  to  produce  an  estimate  of  sensor  positions.  This  offers  several 
advantages,  many  of  which  are  also  offered  by  our  LaSLAT  technique. 
NBP  allows  the  use  of  non-Gaussian  measurement  models.  It  also  pro¬ 
duces  an  uncertainty  measure  that  provides  a  context  for  use  of  its 
sensor  position  estimates. 

A  number  of  authors  [7,  21,  22,  24,  31]  have  proposed  algorithms 
based  on  coordinate  system  stitching.  These  algorithms  follow  a  three 
step  divide-and-conquer  process.  First,  the  network  is  split  into  small 
overlapping  subregions.  Next,  each  subregion  computes  a  local  map , 
which  is  simply  an  embedding  of  the  subregion  in  a  relative  coordinate 
system.  Finally,  adjacent  subregions  are  registered  into  a  common  co¬ 
ordinate  system  using  the  overlapping  nodes  between  local  maps.  By 
performing  registration  recursively,  all  the  subregions  are  incorporated 
into  a  single  global  coordinate  system. 

Moore  et  al.  [22]  suggest  a  particularly  sophisticated  approach  to 
subregion  formation,  in  which  the  subregions  are  chosen  based  on  their 
likelihood  of  forming  accurate  local  maps.  This  technique  has  the  ad¬ 
vantage  of  producing  a  more  homogeneous  accuracy  level  across  the 
network.  However,  it  achieves  this  by  refusing  to  localize  sensors  that 
cannot  be  positioned  using  pair-wise  ranges  to  nearby  nodes. 

This  illustrates  a  fundamental  problem  with  localization  based  on 
inter-sensor  ranges:  the  limited  availability  of  inter-sensor  range  data 
means  that  some  sensors  may  prove  impossible  to  accurately  localize. 
The  advantage  of  simultaneous  localization  and  tracking  (SLAT)  is  that 
these  nodes  may  become  localizable  in  the  future  after  additional  sight¬ 
ings  of  the  mobile.  Furthermore,  SLAT  offers  highest  accuracy  where 
mobiles  move  most  frequently,  whereas  inter-sensor  range-based  algo¬ 
rithms  typically  offer  highest  accuracy  where  sensor  density  is  greatest 
or  in  areas  where  ranging  between  sensors  is  easiest.  In  an  ideal  world, 
these  regions  would  coincide,  but  in  an  ad-hoc  deployment,  this  may 
not  be  the  case. 
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2.2  Tracking 

Like  localization,  target  tracking  is  a  well  established  problem  in  sen¬ 
sor  networks.  However,  most  literature  on  the  subject  assumes  that 
localization  information  is  available  before  tracking  begins.  As  a  re¬ 
sult,  most  of  the  literature  on  tracking  considers  cluster  formation  [9], 
power  economy  [36],  or  related  problems  such  as  target  identification 
and  classification  [3]. 

In  this  thesis,  we  assume  that  targets  are  easy  for  sensors  to  iden¬ 
tify,  allowing  us  to  focus  purely  on  tracking.  This  is  reasonable  in  any 
circumstance  where  the  targets  are  tracked  voluntarily  -  for  instance 
when  the  network  is  guiding  a  mobile  robot  or  providing  in-building 
location  information  to  a  hand-held  computer  user.  We  also  ignore 
concerns  such  as  power  conservation,  and  focus  instead  on  what  in¬ 
formation  can  be  extracted  from  sensor  measurements  once  they  are 
obtained. 


2.3  Calibration 

Several  authors  [6, 35]  have  researched  automatic  calibration  in  sensor 
networks.  Bychkovskiy  et  al.  [6]  calibrate  photo-sensors  in  a  dense 
network  by  making  use  of  the  fact  that  adjacent  sensors  are  likely  to 
observe  the  same  level  of  stimulus. 

The  most  directly  relevant  work  to  this  thesis  is  the  Calamari  Ad- 
hoc  Localization  System  by  Cameron  Whitehouse  [35].  Calamari  uses 
acoustic  time  difference  of  arrival  (TDoA)  ranging  between  sensors  to 
perform  localization.  Whitehouse  showed  that  post-deployment  cali¬ 
bration  of  the  acoustic  sensors  dramatically  improved  localization  ac¬ 
curacy.  He  also  developed  a  macro-calibration  technique  in  which  in¬ 
dividual  sensor  parameters  are  chosen  to  optimize  the  performance  of 
the  network  as  a  whole. 


2.4  SLAT 

Various  authors  have  used  mobiles  to  localize  sensor  networks  [8, 12, 
26-28] ,  but  these  methods  assume  the  location  of  the  mobile  is  known. 
One  exception  is  [12],  which  builds  a  constraint  structure  as  measure¬ 
ments  become  available.  Compared  to  [12],  we  employ  a  very  extensible 
statistical  framework  that  allows  more  realistic  measurement  models. 
Priyantha  et  al.  [28]  show  how  to  guide  a  mobile  to  form  rigidly  local- 
izable  structures.  Our  problem  differs  from  the  one  discussed  in  [28]  in 
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that  we  make  no  assumptions  about  the  trajectory  of  the  mobile. 

Our  method  is  most  similar  to  [26] ,  which  used  an  Extended  Kalman 
Filter  (EKF)  to  track  an  underwater  vehicle  while  localizing  sonar  bea¬ 
cons  capable  of  measuring  their  range  to  the  vehicle.  We  replace  the 
EKF’s  approximate  measurement  model  with  one  based  on  Laplace’s 
method.  This  provides  faster  convergence  and  greater  estimation  ac¬ 
curacy.  We  also  demonstrate  that  the  Bayesian  filtering  framework 
can  calibrate  the  sensor  nodes,  and  that  the  computation  is  capable  of 
distributing  over  the  sensor  nodes  in  a  straightforward  way. 


2.5  SLAM 

The  SLAT  problem  is  similar  to  SLAM  and  the  3D  Structure  from  Mo¬ 
tion  (SFM)  problem  in  computer  vision.  In  these  problems,  two  sets 
of  unknown  variables  are  coupled  in  such  a  way  that  jointly  estimat¬ 
ing  the  two  sets  is  relatively  difficult,  while  estimating  one  set  with 
the  other  set  given  is  relatively  easy.  In  sensor  network  localization, 
knowing  the  position  of  the  nodes  significantly  simplifies  tracking,  and 
knowing  the  position  of  the  mobile  significantly  simplifies  localization. 
In  SLAM  and  SFM,  this  relationship  holds  between  the  pose  of  the 
robot  or  camera  and  that  of  features  in  the  scene. 

Our  solution  to  SLAT  adopts  various  important  refinements  to  the 
original  Extended  Kalman  Filter  (EKF)  formulation  of  SLAM  [33]. 
LaSLAT  processes  measurements  in  small  batches  and  discards  vari¬ 
ables  that  are  no  longer  needed,  as  demonstrated  by  McLauchlan  [19]. 
Following  [34],  LaSLAT  operates  on  inverse  covariances  of  Gaussians 
rather  than  on  covariances  directly  to  accelerate  updates  and  facilitate 
distributed  computation. 
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Chapter  3 


LaSLAT 


LaSLAT  uses  a  Bayesian  filtering  framework.  Under  this  framework, 
as  each  batch  of  measurements  becomes  available,  it  is  used  to  up¬ 
date  a  prior  distribution  over  sensor  locations,  the  mobile  trajectory, 
and  various  sensing  parameters.  The  resulting  posterior  distribution  is 
then  propagated  forward  in  time  using  a  dynamics  model  to  make  it  a 
suitable  prior  for  use  with  the  next  batch  of  measurements. 

In  LaSLAT,  after  incorporating  each  batch  of  measurements  the 
posterior  distribution  is  approximated  with  a  Gaussian  using  Laplace’s 
method  [13].  Consequently,  the  amount  of  state  saved  between  batches 
is  constant  with  respect  to  the  number  of  measurements  taken  in  the 
past.  The  Gaussian  approximation  also  simplifies  propagation  using  the 
dynamics  model  and  incorporation  of  the  next  batch  of  measurements. 


3.1  Approximate  Bayesian  Filtering  for  SLAT 

As  the  mobile  moves  through  the  network,  it  periodically  emits  events 
which  allow  some  of  the  sensors  to  measure  their  distances  from  the 
mobile. 

Let  e3  denote  the  location  of  the  mobile  when  it  generated  the 
jth  event.  The  tth  batch  e*  is  a  collection  of  consecutive  events  e*  = 

{em  . .  ,em+n},  with  e*  denoting  the  jth  event  in  the  fth  batch.  Each 
LaSLAT  iteration  incorporates  the  measurements  from  a  single  batch 
of  events. 

Let  Si  =  [s?  s®]  represent  the  unknown  parameters  of  sensor  i, 

with  denoting  the  sensor’s  position  and  s®  its  calibration  parameters. 
Then  s  =  {s^ }  is  the  set  of  all  sensor  parameters.  The  scalar  j/L  denotes 
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the  range  measurement  between  sensor  i  and  the  jth  event  in  batch  t, 
so  y4  =  {y\j}  is  the  collection  of  all  range  measurements  in  batch  t. 

For  each  batch  t,  e4  and  s  are  the  unknown  values  that  must  be 
estimated.  We  aggregate  these  unknowns  into  a  single  variable  x4  = 
[s  e4]  for  notational  simplicity.  Note  that  as  defined,  e  and  sx  may 
be  vectors  of  either  2D  or  3D  points,  allowing  LaSLAT  to  be  run  easily 
in  either  two  or  three  dimensions. 

The  Bayesian  filtering  framework  is  a  non-linear,  non-Gaussian  gen¬ 
eralization  of  the  Kalman  Filter.  For  each  batch  t,  it  computes  the 
posterior  distribution  over  sensor  parameters,  s,  and  events  locations, 
e4,  taking  into  account  all  range  measurements  taken  so  far: 

P(xV, y2,  •••,/)• 

In  LaSLAT,  we  wish  to  update  this  distribution  as  range  measure¬ 
ments  become  available,  and  discard  measurements  as  soon  as  they 
have  been  incorporated.  To  do  this,  one  can  rewrite  the  distribution  in 
terms  of  a  measurement  model  and  a  prior  distribution  derived  from 
the  results  of  the  previous  iteration.  Rewriting  ^(x^y^y2, . . .  ,y4)  as 
p(x4|yoW,y4),  we  get  by  Bayes’  rule: 

p(xt\yold,yt)<xp(yt,xt\yold) 

ocp(yV,y°WMxVW) 

ocp(yV)p(xVW),  (3.1) 

where  proportionality  is  with  respect  to  x4.  The  final  equality  follows 
because  when  the  sensor  and  mobile  locations  are  known,  the  past  mea¬ 
surements  do  not  provide  any  additional  useful  information  about  the 
new  batch  of  measurements.  The  distribution  p(y4|x4)  is  the  measure¬ 
ment  model:  it  reflects  the  probability  of  a  set  of  observations  given  a 
particular  configuration  of  sensors  and  event  locations  (Section  3.2). 

The  distribution  p(x4  \yold)  summarizes  all  information  collected 
prior  to  the  current  batch  of  measurements,  in  the  form  of  a  predic¬ 
tion  of  x4  and  an  uncertainty  measure.  It  can  be  computed  from  the 
previous  estimate,  p(x4_1  \yold),  by  applying  a  dynamic  model: 

p(x4|yoM)=  [  p(xt-1|y0,<,)p(xt|x4-1)dx4-1.  (3.2) 

Jx*-1 

The  distribution  p(x4|x4-1)  models  the  dynamics  of  the  configuration 
from  one  batch  to  another  by  discarding  old  event  locations  and  pre¬ 
dicting  the  locations  of  new  events  (Section  3.4). 
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1.  Observe  a  new  batch  of  measurements  y*. 

2.  Represent  the  posterior  p(x4| y4,yoid)  in  terms  of  the  prior 
p(x4  |yoW)  and  the  measurement  model  p(y4|x4)  using  Equation 

(3-1). 

3.  Using  Newton-Raphson  [14],  compute  curvature  at  the  mode  of 
p(x4|y4,yoW)  and  use  it  to  construct  the  approximate  posterior 
<7(x4 |y* ,  yoW).  This  posterior  is  the  estimate  for  the  batch  t  (Sec¬ 
tion  3.3). 

4.  Compute  the  prediction  p(x4+1|y4,  yoW)  using  ^(x^y4,  yoW)  (Sec¬ 
tion  3.4). 

5.  Using  the  prediction  as  the  new  prior,  return  to  step  1  to  process 
batch  t+  1. 


Table  3.1:  One  iteration  of  LaSLAT.  Incorporates  batch  t  and  prepares 
to  incorporate  batch  t  +  1. 


When  the  measurement  model  p(y4|x4)  is  not  Gaussian,  the  up¬ 
dates  (3.1)  and  (3.2)  become  difficult  to  compute.  We  handle  the  non- 
Gaussianity  of  the  measurement  model  by  approximating  the  poste¬ 
rior  p(x4 \yold)  with  a  Gaussian  distribution  q(xt \yold)  using  Laplace’s 
method  (Section  3.3).  This  Gaussian  becomes  the  basis  for  the  prior 
distribution  for  the  next  batch,  q  is  much  simpler  to  save  between 
batches  than  the  full  posterior  -  in  particular,  it  allows  all  the  old  mea¬ 
surements  to  be  discarded.  Table  3.1  summarizes  the  steps  of  LaSLAT. 

Other  approximate  Bayesian  filters  such  as  the  Extended  Kalman 
Filter  (EKF)  or  particle  filters  could  also  be  used  in  place  of  our  Lapla- 
cian  method.  The  EKF  differs  from  our  algorithm  because  it  does  not 
perform  a  full  optimization  when  incorporating  each  event.  In  many 
cases  this  is  a  helpful  optimization;  however,  as  we  show  in  chapter  6, 
on  the  SLAT  problem  it  sacrifices  accuracy  and  speed  of  convergence. 
Particle  filters  allow  a  closer  approximation  of  the  posterior  distribu¬ 
tion,  especially  when  the  distribution  is  multi-modal.  However,  our  al¬ 
gorithm  seems  to  perform  well  in  practice  while  requiring  significantly 
less  computation. 


18 


3.2  Measurement  Model 


New  measurements  influence  localization  and  calibration  estimates  via 
the  measurement  model.  A  measurement  model  is  a  probability  distri¬ 
bution  ,p(yt  x*)  over  a  batch  of  range  measurements,  given  a  particular 
choice  of  the  calibration  parameters  and  positions  for  the  sensors  and 
the  mobile.  One  advantage  of  the  LaSLAT  framework  is  that  this  mea¬ 
surement  model  can  be  tailored  to  specific  measurement  hardware  and 
deployment  parameters.  The  measurement  model  encapsulates  details 
such  as  the  kind  of  data  being  observed  (angle  of  arrival,  distance,  radio 
signal  strength,  etc),  as  well  as  the  types  of  noise  that  are  possible  in 
the  environment. 

In  this  thesis,  we  develop  a  measurement  model  that  reflects  the 
characteristics  of  the  popular  acoustic  time  difference  of  arrival  (TDoA) 
ranging  technique.  In  most  TDoA  implementations,  the  transmitter  (in 
this  case  the  mobile)  emits  a  tagged  radio  message,  then  a  short  time 
later  produces  an  acoustic  signal.  Upon  hearing  the  radio  message, 
the  sensors  in  the  network  activate  their  microphones  and  listen  for 
the  arrival  of  the  acoustic  signal.  Once  the  acoustic  signal  arrives,  the 
sensors  use  the  difference  in  arrival  time  between  the  acoustic  signal 
and  the  radio  message  to  estimate  the  distance  between  the  sensor  and 
the  transmitter. 

This  technique  can  be  quite  accurate,  especially  when  line-of-sight 
exists  between  the  transmitter  and  the  receiver.  The  Cricket  ranging 
system  [1],  which  uses  a  single  inaudible  ultrasound  pulse,  can  measure 
distances  up  to  ten  meters  with  error  as  small  as  1-2  centimeters. 

However,  TDoA  measurements  are  susceptible  to  several  types  of 
error.  First,  various  random  delays  in  the  ranging  process  introduce 
small  (sub-centimeter)  errors,  which  take  the  form  of  a  variance  from 
measurement  to  measurement.  Second,  mis-calibration  can  result  in  a 
measurement  bias.  For  instance,  temperature  or  humidity  can  cause 
the  actual  speed  of  sound  to  be  different  from  the  pre-calibrated  con¬ 
stant  used  by  the  sensor.  The  measurement  bias  can  easily  vary  from 
deployment  to  deployment,  making  accurate  pre-calibration  at  the  fac¬ 
tory  nearly  impossible.  In  the  Cricket  ranging  system,  we  have  ob¬ 
served  biases  of  between  3  and  10  centimeters.  The  biases  vary  only 
slightly  from  sensor  to  sensor,  but  are  typically  consistent  throughout 
a  deployment  area.  Finally,  TDoA  measurements  can  be  vulnerable  to 
both  echoes  and  ambient  noise.  In  both  cases  errors  occur  because  the 
time  of  flight  observed  by  the  sensor  does  not  correspond  to  the  straight 
line  distance  from  the  sensor  to  the  event;  instead,  it  corresponds  to  a 
longer  path  (in  the  case  of  an  echo),  or  no  path  at  all  (in  the  case  of  an 


19 


ambient  noise).  In  one  test  environment,  a  security  system  produced  40 
kHz  ultrasound  that  occasionally  interfered  with  the  Crickets’  ranging 
system. 

In  this  section,  we  ignore  echo  effects,  and  assume  that  each  mea¬ 
surement  is  a  corrupted  version  of  the  true  distance  between  the  event 
and  the  sensor  that  took  the  measurement: 

y\j  =  llsi  ~  ejll  +  si  +  u\j,  (3.3) 

where  |  •  ||  indicates  the  vector  2- norm,  giving  the  Euclidean  distance 
between  s?  and  e( .  is  a  zero- mean  Gaussian  random  variable  with 

1  J  LJ 

variance  er2,  and  s®  is  a  bias  parameter  that  models  an  unknown  shift 
due  to  sensor  mis-calibration.  In  section  4.1  we  present  a  richer  model 
that  performs  better  in  the  presence  of  echoes. 

As  defined  in  equation  (3.3),  p(yL|s,,  e*)  is  a  univariate  Gaussian 
with  mean  ||si  —  e*-||  +  sf  and  variance  er2.  Since  each  measurement 
y*.  depends  only  on  the  sensor  s,  that  took  the  measurement  and  the 
location  e*  of  the  mobile  when  it  generated  the  event,  the  measurement 
model  factorizes  according  to 

piy*^)  =  Y[p{ytij\si,etj),  (3.4) 

i,j 

where  the  product  is  over  the  sensors  and  the  events  that  they  perceived 
in  batch  t.  Equation  (3.4)  is  the  complete  measurement  model  for  a 
batch  of  measurements. 

Though  LaSLAT  is  designed  to  operate  in  their  absence,  LaSLAT 
can  make  use  of  inter-sensor  range  measurements  if  they  are  available. 
These  can  be  encoded  using  a  generative  model  similar  to  equation 
(3.3): 

Vkl  =  llsfc  —  s?  II  +  sk  +  (3-5) 

where  s^.  and  s /  are  the  sensor  pair  that  generates  the  inter-sensor  range 
measurement  y*kl,  and  is  additive  zero-mean  Gaussian  noise.  The 
resulting  probability  distribution  p(ytk[\sk,si)  becomes  an  additional 
factor  in  the  product  (3.4).  This  demonstrates  the  ease  with  which 
LaSLAT  can  rigorously  incorporate  additional  sources  of  localization 
information. 


3.3  Incorporating  Measurements 

The  measurement  model  derived  in  the  last  section  can  be  combined 
with  a  prior  distribution  ^(x^y0^)  using  equation  (3.1)  to  find  the 
posterior  distribution  p(xt\yold ,yt). 
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This  posterior  has  no  compact  representation  -  in  particular,  it 
takes  up  space  proportional  to  the  total  number  of  measurements  ob¬ 
served  by  the  network  since  the  first  batch.  To  curb  this  complexity, 
we  save  only  a  Gaussian  approximation  of  the  posterior.  This  Gaussian 
representation  requires  constant  space  with  respect  to  the  number  of 
measurements  observed,  which  allows  LaSLAT  to  run  indefinitely  with¬ 
out  requiring  additional  memory.  It  also  has  the  advantage  of  being 
easy  to  distribute  across  the  sensor  network  (chapter  5). 

This  approximate  Gaussian  posterior  g(x4  |y°  ,  y4)  can  be  obtained 
from  the  prior  distribution  p(x4 \yold)  and  the  measurement  model  p(yt  |x4) 
using  Laplace’s  method  [13]. 

To  fit  an  approximate  Gaussian  distribution  q(x )  to  a  distribution 
p(x),  Laplace’s  method  first  finds  the  mode  x*  of  p(x),  then  computes 
the  curvature  of  the  negative  log  posterior  at  x*. 

Q 2 

A  1  =  -^2l0S^)x_x,- 

The  mean  and  covariance  of  q(x)  are  then  set  to  x*  and  A  respectively. 
Notice  that  when  p  is  Gaussian,  the  resulting  approximation  q  is  exactly 
p.  For  other  distributions,  the  Gaussian  q  locally  matches  the  behavior 
of  p  about  its  mode. 

The  mode  finding  problem  can  be  expressed  as: 
x4*  =  argmaxp(x4|y4,  yold) 

X* 

=  arg  min  -  log  [p(yt|xt)p(xt|yoW)] 

=  arg  min  (x4  -  /r)Tffx(x4  -  p)  +  ^  JI(||sf  -  e4 1|  +  s  ■  -  y^)2, 

(3.6) 

where  p  =  E  [x4|y4,yoW] ,  and  Hx  =  Cov-1  [x4|yoW] . 

We  use  the  Newton-Raphson  iterative  optimization  algorithm  [14] 
to  find  the  mode  x4*  and  the  curvature  H  (Appendix  B).  Following 
Laplace’s  method,  the  mean  E  [x4|yoW,  y4]  of  q  is  set  to  x4*  and  its  in¬ 
verse  covariance  Cov-1  [x4|yoid,y4]  is  set  to  H.  Representing  q  using 
its  inverse  covariance  allows  us  to  avoid  computing  the  matrix  inverse 
H-1  after  adding  each  measurement,  which  significantly  improves  per¬ 
formance  and  facilitates  a  distributed  implementation  of  our  algorithm 
(Chapter  5). 

The  Gaussian  approximation  described  in  this  section  works  well 
when  the  posterior  distribution  has  a  single  strong  mode.  However, 
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when  the  posterior  contains  several  equally  probable  modes,  the  Gaus¬ 
sian  approximation  effectively  discards  all  but  one.  Consequently,  this 
algorithm  is  vulnerable  to  substantial  errors  when  the  measurements 
and  prior  do  not  favor  a  unique  estimate.  This  can  occur  if  the  mobile 
does  not  move  very  much  during  a  batch  and  the  prior  distribution  is 
mostly  uninformative.  In  that  case,  the  Gaussian  approximation  risks 
committing  too  quickly  to  a  particular  parameter  estimate.  In  our  ex¬ 
perience,  this  problem  can  be  avoided  by  using  larger  batch  sizes  when 
the  prior’s  covariance  is  large.  It  can  also  be  solved  by  using  a  Gaussian 
mixture  model  or  a  particle  filter  to  approximate  the  posterior  distribu¬ 
tion  between  batches.  These  alternatives  are  able  to  fit  a  multi-modal 
posterior;  however,  they  are  more  complex  and  consequently  require 
greater  computational  power  to  manipulate. 


3.4  Dynamics  Model 

In  this  thesis,  we  assume  mobiles  can  move  arbitrarily  and  that  sensors 
are  stationary.  When  propagating  the  posterior  g(x4 |yoid,y4)  forward 
in  time,  we  need  only  retain  the  components  that  are  useful  for  incor¬ 
porating  the  next  batch  of  measurements.  Thus,  we  may  remove  the 
estimate  of  the  mobile’s  trajectory  from  batch  t.  but  we  must  incorpo¬ 
rate  a  guess  for  the  mobile’s  path  during  batch  t  +  1.  Therefore,  the 
prediction  step  of  Equation  (3.2)  can  be  written: 

p(xt+1|yoW,yt)  =p(s,e4+Vw,y4)  =  p(et+1)q(s\yMy)  (3.7) 

<?(s|yoW,yt)=  [  g(x4|yoW,y4)<ie4. 

Je * 

The  Gaussian  q(xt\yold,  y4)  captures  the  posterior  distribution  over  sen¬ 
sor  locations  given  all  measurements  taken  so  far,  and  has  already  been 
computed  by  the  method  of  section  3.3.  We  obtain  q(s\yold,  y4),  by 
marginalizing  out  the  mobile’s  trajectory  during  batch  t. 

The  prior  p(e4+1)  is  Gaussian  with  very  broad  covariance,  indicating 
that  the  future  trajectory  of  the  mobile  is  unknown.  In  some  applica¬ 
tions,  it  may  be  possible  to  use  past  trajectories  to  make  better  guesses 
for  e4+1.  For  example,  p(e4+1)  could  be  used  to  require  events  to  form  a 
smooth  path.  This  can  help  accurately  position  events  with  few  quality 
measurements.  For  maximum  generality,  we  will  not  attempt  to  do  so 
in  this  section,  meaning  that  the  mobile  is  allowed  to  move  arbitrarily 
between  events.  In  section  4.2,  we  present  a  method  of  requiring  that 
events  follow  a  smooth  trajectory. 
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The  operations  of  Equation  (3.7)  can  be  carried  out  numerically 
by  operating  on  the  mean  and  inverse  covariance  of  (7(x4|yoid).  First, 
partition  according  to  s  and  e4: 


E  [x4|yoW 


'  £[S|y°id] 
E[e* \yold] 


Cov”1  [x4|yokr 


-  Is  Use4 

Hgt  g  Hgt 


Marginalizing  out  e4  produces  a  distribution  g(s|yoid)  whose  mean  is 
the  s  component  of  the  mean  of  q(xt\yold)  and  whose  inverse  covariance 
is: 

Cov-1  [s|yoW]  =  fis  -  fise^e^^s-  (3-8) 

The  parameters  of  p(x4+1  |yoW)  are  those  of  q( s|yoW),  augmented  by 
zeros  to  account  for  an  uninformative  prior  on  e4+1: 


E  [x4+1|yoM] 
Cov"1  [x4+1|yoW] 


E  [s|yoW]" 

0 

Cov-1  [s|yoW]  0 

0  0 


(3.9) 


The  components  of  the  inverse  covariance  of  p(x4+1  \yold)  corresponding 
to  e4+1  are  set  to  0,  corresponding  to  infinite  variance,  which  in  turn 
captures  our  lack  of  a  priori  knowledge  about  the  location  of  the  mobile 
in  the  new  batch.  The  mean  is  arbitrarily  set  to  0.  If  some  information 
is  known  a  priori  about  e4+1,  then  the  0  components  of  E  [x4+1|yoid] 
and  the  bottom  right  0  components  of  Cov-1  [x4+1|yoW]  can  be  used 
to  capture  that  knowledge. 


3.5  Prior  Information  and  Initialization 

Prior  information  about  the  sensor  parameters  is  easy  to  incorporate 
into  LaSLAT.  Such  information  might  be  available  because  the  sensors 
were  placed  in  roughly  known  positions,  or  because  another  less  accu¬ 
rate  source  of  localization  is  available.  Prior  information  may  also  be 
used  to  define  LaSLAT’s  global  coordinate  system,  by  fixing  the  rela¬ 
tive  positions  of  several  sensors.  In  addition,  calibration  in  the  factory 
might  supply  prior  information. 

If  such  prior  information  is  available  it  can  be  supplied  as  the  prior 
when  incorporating  the  first  batch  of  measurements.  We  set  the  co- 
variance  of  this  prior  to  (Jof  with  cq  a  large  scalar,  which  makes  the 
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prior  diffuse.  The  large  covariance  allows  measurements  to  override 
the  positions  prescribed  by  the  prior,  but  provides  a  sensible  default 
when  few  measurements  are  available.  The  mode  of  this  prior  (or  for 
subsequent  iterations,  the  mode  of  p(s\yold ))  is  also  used  as  the  initial 
iterate  for  the  Newton-Raphson  iterations.  To  obtain  the  initial  iterate 
for  an  event,  we  use  the  average  estimated  location  of  the  three  sensors 
with  the  smallest  range  measurements  to  the  event. 

In  our  experiments,  we  utilize  the  radio  connectivity  of  the  sen¬ 
sors  to  obtain  prior  localization  information.  The  initialization  step 
described  by  Priyantha  et  al.  [29]  (see  Appendix  A)  provides  rough 
position  estimates  to  serve  as  a  prior  before  any  measurements  are  in¬ 
troduced.  This  prior  takes  the  formp(sa:)  <x  exp  ||sf  —  a?°||2  , 

where  x ?  is  the  position  of  the  ith  sensor  as  predicted  by  the  initializa¬ 
tion  step  and  00  is  a  large  variance. 

We  have  observed  empirically  on  Cricket  hardware  [1]  that  though 
measurement  biases  vary  between  deployments,  they  tend  to  be  fairly 
consistent  from  sensor  to  sensor.  This  information  can  be  used  as  a 
prior  over  the  sensor  measurement  biases  sf .  We  encode  this  informa¬ 
tion  as  a  distribution  with  the  form  p(se)  oc  exp  —  J2i~j  isi  —  sj)2  > 
where  the  summation  is  over  sensors  that  are  in  close  proximity  to 
each  other,  and  a &  is  a  constant  used  to  tune  the  impact  of  this  prior 
on  the  measurement  biases.  This  illustrates  a  important  attribute  of 
LaSLAT:  it  is  easy  to  add  platform-specific  optimizations  to  improve 
performance.  In  many  localization  algorithms  this  is  impossible:  a 
change  in  hardware  or  deployment  environment  often  necessitates  a 
completely  new  approach. 
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Chapter  4 

LaSLAT  Extensions 


In  the  last  chapter,  we  described  the  core  LaSLAT  algorithm.  In  this 
chapter,  we  consider  several  useful  enhancements.  First,  we  describe 
a  measurement  model  that  provides  robust  detection  and  elimination 
of  erroneous  range  measurements.  Then,  we  show  how  to  require  that 
LaSLAT’s  event  position  estimates  follow  a  smooth  path. 


4.1  Robust  Measurement  Outlier  Rejection 


When  using  time  difference  of  arrival  (TDoA)  acoustic  ranging,  many 
external  factors  can  cause  non-Gaussian  measurement  error.  For  in¬ 
stance,  in  some  environments,  ambient  sounds  can  cause  completely 
spurious  measurements.  In  addition,  environmental  echoes  can  cause 
substantial  delays  in  the  sound  pulse’s  time  of  flight,  resulting  in  highly 
errant  measurements.  This  effect  is  illustrated  in  figure  4.1.  In  these 
types  of  environments,  it  can  be  helpful  to  define  a  more  robust  mea¬ 
surement  model.  We  propose  a  mixture  model  for  this  purpose.  Let 
h\j  be  a  Bernoulli  random  variable  that  is  1  with  probability  p.  Then 
a  generative  model  of  measurements  might  be: 


—  e" 


4  =  1 
4  =  o 


(4.1) 


where  u\j  is  a  uniform  random  variable  over  all  possible  range  mea¬ 
surements.  Equation  (4.1)  produces  an  accurate  measurement  with 
probability  p1  and  an  uninformative  and  therefore  useless  result  with 
probability  1  —  p.  Physically,  when  a  measurement  occurs,  there  are 
three  possible  results: 
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Reported  distance  (cm) 


Figure  4.1:  An  example  of  outlying  range  measurements.  This  data 
was  generated  by  rotating  a  beacon  cricket  at  a  fixed  distance  from 
a  sensor.  This  histogram  plots  the  range  measurements  taken  by  the 
sensor.  Notice  that  the  majority  of  the  measurements  are  reasonably 
close  to  the  true  distance  (much  of  the  variance  is  caused  by  the  move¬ 
ment  of  the  ultrasound  transmitter  on  the  beacon  with  respect  to  the 
sensor) ;  however,  a  few  measurements  are  nearly  half  a  meter  too  long. 
These  outlying  measurements  can  dramatically  decrease  localization 
and  tracking  performance,  so  it  is  desirable  to  detect  and  discard  them. 
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•  The  radio  message  is  received,  and  the  sensor’s  microphone  is  ac¬ 
tivated.  However,  before  the  acoustic  pulse  arrives  at  the  sensor, 
an  unexpected  environmental  sound  reaches  the  sensor,  causing 
a  short  measurement. 

•  The  radio  message  is  received,  and  the  sensor’s  microphone  is 
activated.  However,  the  acoustic  pulse  is  deflected  by  an  obstacle. 
Before  the  microphone  deactivates  due  to  a  timeout,  an  echo  of 
the  acoustic  pulse  arrives,  causing  a  long  measurement. 

•  The  radio  message  is  received,  and  the  acoustic  pulse  arrives  at 
the  microphone  without  incident.  The  measurement  is  accurate. 


Equation  (4.1)  models  this  as  a  Bernoulli  trial:  the  third  case  occurs 
with  probability  p ,  and  one  of  the  error  cases  occurs  with  probability 
1  —  p.  h\j  is  the  random  variable  representing  the  outcome  of  this  trial 
for  measurement  y\j.  Equation  (4.1)  leads  directly  to  the  probability 
distribution: 


PiUij  I  Si  ,  )  p(/r,;.y 


1 )  vy-  exp 


(ulj  —  IIs?  —  ej  II  —  si)2 

2a2 


+  V{h\j 


(4.2) 


where  Z\  and  Zi  are  normalization  constants.  This  changes  the  struc¬ 
ture  of  equation  (3.6),  since  each  measurement  y|  -  now  has  an  unknown 
explanatory  variable  hy.  This  new  structure  is  straightforward  to  op¬ 
timize  using  an  Expectation  Maximization  (EM)  algorithm  [13].  Let 
9  =  [x*].  EM  finds  a  maximizing  9  repeatedly  applying  the  following 
update  step: 


9k+1  =  argmax  Ep{ht [logp(y\ h'|0, yold)  +  \ogp(9\yold)\ 

u 

=  argmax  ^p(h;|y*,  9k)  [logp(y‘,  h* |6»,  yold)  +  \ogp(9\yold)] 

h* 

=  argmax  \ogp(9\yold)  +  EE  p(hij\ytj,  ok)  !og  PivljMj \o) 

(4.3) 

4.1.1  E-step 

In  the  E-step,  we  must  compute  p{hy  \y\p  6k)  for  all  i  and  j.  Since 
log p{y\j,hy  =  O|0)  is  constant  with  respect  to  9,  the  corresponding 
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terms  have  no  effect  on  the  maximization  (4.3).  Consequently,  we  need 
only  calculate: 


P(h\j 


i|  ytj,ok) 


P{jjIpK3  =  ijgfc) 
'EhP(Vij>h\ok) 

P(htij  =  l)p(ytij\hjj  =  l,0k) 
E  hv{h)v{ylj\h,ok) 


This  expression  is  easy  to  compute  since  it  is  in  terms  of  the  measure¬ 
ment  model  conditioned  on  h  and  the  Bernoulli  distribution  p(h): 


V{h\j 


pN 

pN  +  (1  -p)/Z2' 


(4.4) 


where 


Mi 


IIs?1  ~  ej II  ~  si)2 

2  a2 


4.1.2  M-step 

In  the  M-step,  we  optimize  equation  (4.3)  to  find  a  new  estimate  9k+1 . 
Let  w\'*  =  p{h\j  =  1| ylj,9k)  as  computed  in  the  E-step  by  eciuation 
(4.4).  Equation  (4.3)  reduces  to: 

9k+1  =  arg max  logp(%oW)  +  ^  w\f  logp(yb,  h\j  =  1\0) 

ij 

=  arg  min  (x*  —  /r)Tnx(xt  —  p) 

+  llSi'  “  e‘H  +  Si  -  ylif-  (4-5) 

The  M-step  (4.5)  is  therefore  a  re-weighting  of  the  original  LaSLAT 
optimization  problem  (equation  (3.6)).  As  a  result,  it  can  be  per¬ 
formed  analogously  using  Newton- Raplrson  (see  Appendix  B).  Note 
that  the  weight  tub  of  a  measurement  t/b  after  the  final  EM  iteration 
corresponds  to  the  probability  that  the  measurement  is  accurate  given 
the  current  parameter  estimates.  Thus,  accurate  measurements  are  as¬ 
signed  weights  close  to  one,  and  highly  inaccurate  measurements  are 
assigned  weights  close  to  zero.  This  accomplishes  the  goal  of  rejecting 
outlying  measurements. 
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4.1.3  Outlier  rejection  summary 

When  the  update  step  (4.3)  is  performed  repeatedly,  the  9k's  con¬ 
verge  to  x4*,  an  optimal  estimate  of  sensor  parameters  and  event  lo¬ 
cations  that  detects  and  ignores  outlying  measurements.  This  im¬ 
proves  LaSLAT’s  performance  in  the  presence  of  environmental  am¬ 
bient  sounds  or  echoes  due  to  physical  obstacles. 

4.2  Specifying  Mobile  Dynamics 

In  section  3.4,  we  assumed  that  the  mobile  moved  arbitrarily.  However, 
in  some  cases  (for  instance  when  tracking  targets  that  move  continu¬ 
ously  and  trigger  events  frequently),  one  may  improve  performance  by 
explicitly  requiring  that  successive  events  occur  close  to  each  other.  In 
this  section,  we  present  a  technique  for  adding  such  a  constraint  to 
LaSLAT. 

As  shown  in  section  3.4,  the  prior  distribution  for  LaSLAT 
p(x4+1|y4,  yoW)  is  computed  according  to: 

p(xt+1\yoldy)  =p(et+1)q(s\yMy).  (4.6) 

Until  now,  we  have  defined  p(e4+1)  to  be  a  Gaussian  with  high  covari¬ 
ance.  Let  e4+1  =  [e*  ej  e“],  where  ev-  and  e“  represent  the  velocity 
and  acceleration  of  the  mobile  at  the  time  of  event  j,  and  define  eg+1  to 
be  the  estimated  parameters  of  the  event  immediately  preceding  batch 
t  +  1.  Then  we  can  express  p(et+1)  as  follows: 


p(et+1)  =  n^ei+1|ei-1i)- 

j= 1 


We  enforce  smoothness  by  setting  p(e4+1|e4+1)  to  a  Gaussian  with 
mean  and  covariance  cr^I.  A  is  a  matrix  that  expresses  the 

expected  physical  motion  of  the  mobile.  We  assume  that  events  occur 
at  a  constant  rate,  and  ignore  quadratic  terms,  allowing  us  to  use  the 
matrix: 


A  = 


II  0 
0  II 
0  0  1 


as  the  dynamics  matrix.  This  matrix  has  the  advantage  of  being  compu¬ 
tationally  straightforward  while  having  the  desired  effect  of  smoothing 
the  event  positions,  tjd  is  a  tunable  constant  that  allows  us  to  vary  the 
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effect  of  the  smoothness  prior  versus  the  effect  of  measurement  data  on 
the  events’  position  estimates. 

Thus,  we  can  rewrite  equation  (4.7)  as: 


p(et+  )  a  exp 


This  term  can  be  rewritten  in  the  form: 

1 


p(et+  )  a  exp 


~  Pdf'  B{et  -  pd) 


where  the  block  tridiagonal  matrix  B  and  vector  pd  are  constants  ex¬ 
pressible  in  terms  of  A  and  e^+  : 


I  +  ATA  - AT 
-A 

'  •  •  I  +  AT A 
-A 


Pd  =  B  1 


-A 

0 


e 


t+ 1 

o 


-AT 

I 


B  and  vector  pd  can  then  incorporated  into  equation  (3.9)  as  follows: 


E  [xt+1  \yold]  = 

'E  [s  yoW]" 

Pd 

~ 1  [xt+1  \yold]  = 

Cov^1  [s|y 

0 

B 


As  we  demonstrate  in  chapter  6,  this  smoothness  constraint  can 
noticeably  improve  performance  on  occasional  poorly  measured  events. 

Unfortunately,  the  smoothness  constraint  complicates  the  marginal¬ 
ization  step  shown  in  equation  (3.8),  since  the  inverted  matrix  Q is  no 
longer  diagonal.  Consequently,  the  smoothness  prior  is  recommended 
for  use  only  when  LaSLAT  computations  are  being  performed  centrally. 
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Chapter  5 

Implementation 


Our  current  implementation  sends  measurement  batches  to  a  central 
computer;  however,  we  show  here  that  LaSLAT  can  be  feasibly  dis¬ 
tributed  if  desired.  We  also  explore  running  time  bounds  for  both  a 
distributed  implementation  and  a  centralized  implementation.  Finally, 
we  discuss  the  trade-off  between  centralized  and  distributed  processing 
in  LaSLAT. 


5.1  Definitions 

In  order  to  quantify  LaSLAT’s  performance,  several  definitions  are  re¬ 
quired.  Let  niocai  be  the  expected  number  of  sensors  within  a  one-hop 
neighborhood  of  a  sensor  s.  The  one-hop  neighborhood  can  be  thought 
of  as  a  circle  around  s  whose  radius  is  two  times  the  sensor’s  maximum 
sensing  range.  Thus,  the  one-hop  neighborhood  contains  all  sensors 
that  can  sense  an  event  in  common  with  s.  We  assume  that  s  can  com¬ 
municate  directly  with  all  of  the  sensors  in  this  one-hop  neighborhood. 

Furthermore,  we  define  nevents  to  be  the  number  of  events  observed 
by  s  in  a  single  batch  of  measurements.  nevents  can  be  held  constant 
by  varying  the  amount  of  time  between  batches. 

Note  that  niocai  and  nevents  are  constants  fixed  at  deployment  by 
the  network  designer.  This  means  that  asymptotic  bounds  in  niocai  and 
nevents  are  in  some  sense  constant  bounds,  since  they  enable  sensors  to 
be  provisioned  with  an  amount  of  memory  and  processing  power  that  is 
guaranteed  sufficient  no  matter  how  many  events  the  network  observes. 
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5.2  Graph  Locality 

The  Gaussian  prior  p(xt\y0  <  )  is  completely  summarized  by  a  vector  of 
means  x4*  and  an  inverse  covariance  Cl. 

The  symmetric  inverse  covariance  matrix  Cl  =  [  qs  e  ]  defines  an 
undirected  graph  between  sensors  and  events.  Two  vertices  in  this 
graph  are  connected  if  their  corresponding  block  in  Cl  is  non-zero.  We 
say  Cl  has  local  connectivity  if  the  corresponding  graph  only  connects 
sensors  that  are  within  one  hop  of  each  other  and  connects  events  only 
to  the  sensors  that  measured  the  event.  As  we  will  show,  in  LaSLAT 
Cl  always  has  local  connectivity. 


5.3  Distributability  of  the  Prior 

The  mean  vector  is  straightforward  to  distribute:  each  sensor  simply 
stores  its  own  mean,  those  of  its  one  hop  neighbors,  and  those  of  all 
nearby  events.  Thus,  storage  for  means  requires  at  most  0(niocai  + 
Uevents )  per  sensor. 

Cl  requires  more  careful  consideration,  but  is  also  distributable.  If 
Cl  has  local  connectivity,  each  row  of  Cl  corresponding  to  a  sensor  has 
about  niocai  +  nevents  non-zero  entries.  Each  row  corresponding  to  an 
event  has  less  than  niOCai  non-zero  elements,  since  only  sensors  within 
the  neighborhood  of  the  event  obtain  measurements  to  it.  Locally 
connected  matrices  are  therefore  easy  to  distribute.  Each  sensor  stores 
its  own  rows  in  Cl.  Event  rows  are  delegated  randomly  to  a  sensor 
for  storage,  leaving  sensor  storing  about  nevents / niocai  event  rows.  The 
amount  of  data  stored  by  each  sensor  is  consequently  0(niocai+nevents)  ■ 
If  the  computation  is  performed  centrally,  then  the  central  computer 
must  store  this  same  amount  of  data  per  sensor  in  the  network. 


5.4  Performing  Computations 

LaSLAT  consists  of  two  significant  computational  steps:  incorporating 
measurements  and  applying  the  dynamics  model.  As  formulated  in  this 
paper,  these  steps  can  be  performed  using  only  local  communication 
between  sensors  that  have  witnessed  a  common  event.  Furthermore, 
these  operations  retain  local  connectivity  in  the  prior  inverse  covariance 
matrix  Cl. 
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5.4.1  Measurement  incorporation 

The  principal  operation  involved  in  incorporating  new  measurements  is 
a  Newton- Raphson  iteration,  shown  in  equation  (B.3).  Each  Newton- 
Raplison  iteration  has  two  parts.  First  a  matrix  and  vector  must  be 
computed  based  on  equation  (B.3).  Then,  a  least  squares  optimization 
must  be  performed. 

The  matrix  and  vector  can  be  computed  locally,  since  they  require 
only  the  parts  of  fi  and  x**  that  are  found  locally  and  the  measurements 
to  any  local  events.  The  communication  and  time  costs  for  each  are 
proportional  to  nevents  *  niocai ,  which  is  the  minimum  time  required  for 
each  sensor  to  broadcast  new  measurements  to  neighboring  sensors  and 
receive  their  measurements  in  return.  It  is  similarly  the  minimum  time 
for  all  sensors  to  transmit  their  measurements  to  a  central  computer. 

Once  the  matrix  and  vector  are  computed,  the  least  squares  opti¬ 
mization  can  be  performed  using  Gauss-Seidel  iterations  [2],  Gauss- 
Seidel  is  guaranteed  to  converge  when  solving  symmetric  positive  defi¬ 
nite  systems  of  equations  like  those  found  in  LaSLAT.  Each  iteration  of 
Gauss-Seidel  requires  O(niocai)  computation  and  0(1)  radio  messages 
per  sensor  and  event  when  distributed,  or  0{n  *  niocai)  computation 
when  centralized,  where  n  is  the  total  number  of  sensors  and  events. 
In  practice,  we  find  that  Gauss-Seidel  converges  in  a  few  tens  of  iter¬ 
ations  for  our  systems,  since  LaSLAT  does  not  require  high  precision 
convergence. 

Gauss-Seidel  requires  that  sensors  perform  their  processing  in  a  con¬ 
sistent  order,  which  diminishes  the  potential  parallelization  of  the  least 
squares  computation.  However,  with  a  constant  bound  on  the  number 
of  Gauss-Seidel  iterations,  the  total  time  required  for  each  distributed 
Newton-Raphson  iteration  is  O(n*niocai).  Note  that  this  is  not  a  tight 
upper  bound:  as  the  sensor  parameters  begin  to  converge,  many  pa¬ 
rameters  will  not  need  to  be  updated  every  batch.  This  increases  the 
amount  of  parallelism  that  can  be  exploited,  allowing  the  total  running 
time  to  approach  O(n;oco;),  the  amount  of  time  required  to  simply 
locate  the  newest  events.  See  [2]  for  more  details  and  an  in-depth  de¬ 
scription  of  Gauss-Seidel  iterations. 

5.4.2  Event  marginalization 

Marginalization  is  performed  using  equation  (3.8): 

Gov-1  [s|yoid]  =  -  f 'lsetCl~t1Clets. 

It  is  distributable  because  each  sensor  row  is  updated  only  on  behalf  of 
local  events.  Since  fise*  and  flets  are  sparse  and  Het  is  block  diagonal, 
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the  total  time  required  is  only  0(niocai  *  nevents)  per  sensor.  All  the 
computations  can  occur  in  parallel. 

Unfortunately,  the  smooth  dynamics  extension  developed  in  section 
4.2  produces  an  Clet  that  is  block  tridiagonal  and  has  a  dense  matrix 
inverse  in  general.  As  a  result,  the  smooth  dynamics  extension  ruins 
local  connectivity.  Thus,  smooth  dynamics  may  only  be  employed  when 
LaSLAT  computations  are  performed  at  a  central  computer. 

5.4.3  Preservation  of  local  connectivity 

It  remains  to  be  shown  that  the  LaSLAT  computations  retain  local  con¬ 
nectivity  in  Cl  when  the  smooth  dynamics  extension  is  not  employed. 
This  is  easily  confirmed  by  induction  on  Cls.  The  initial  prior  has 
Cls  =  ctqI,  which  is  diagonal  and  therefore  locally  connected.  As  we 
show  in  the  appendix,  the  measurement  incorporation  and  Gaussian 
approximation  steps  do  not  change  the  connectivity  of  Cl.  The  con¬ 
nectivity  of  Cls  only  changes  when  events  are  marginalized  out  of  the 
Gaussian  prediction  by  equation  (3.8).  It  can  be  verified,  however,  that 
the  —  term  added  to  Cls  only  affects  elements  of  fls  whose 

corresponding  sensors  observed  an  event  in  common  during  the  most 
recent  batch.  As  a  result,  Cls  retains  local  connectivity  during  LaSLAT 
operations. 

5.4.4  LaSLAT  convergence 

Each  LaSLAT  batch  requires  at  least  one  Newton-Raphson  optimiza¬ 
tion,  which  consists  of  several  iterations.  However,  these  iterations 
need  not  continue  until  the  solution  is  fully  converged.  In  fact,  if  each 
LaSLAT  batch  performs  only  one  Newton-Raphson  iteration,  then  the 
resulting  algorithm  is  almost  precisely  the  Extended  Kalman  Filter 
(EKF)  form  of  SLAT.  As  we  show  in  chapter  6,  additional  Newton- 
Raphson  iterations  substantially  improve  performance.  However,  little 
performance  is  lost  if  the  number  of  iterations  is  bounded  at  a  small 
constant.  In  our  experiments,  Newton-Raphson  often  converged  in  less 
than  ten  iterations.  In  two  dimensions,  frequently  as  few  as  three  or 
four  were  required  for  convergence.  Thus,  the  number  of  iterations  may 
be  treated  as  a  constant  factor. 

When  using  measurement  outlier  rejection  (section  4.1),  even  less 
Newton-Raphson  iterations  are  required  per  optimization,  since  the  EM 
optimization  runs  Newton-Raphson  repeatedly.  The  EM  optimization 
remains  distributable,  since  the  only  additional  processing  step  is  the 
computation  of  measurement  weights. 
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5.5  Centralized  vs.  Distributed  Implemen¬ 
tation 

As  we  have  demonstrated,  LaSLAT  is  amenable  to  both  centralized 
and  distributed  implementation.  In  a  centralized  implementation,  the 
network  must  transport  all  range  observations  to  a  central  computer, 
which  performs  LaSLAT  computations  and  optionally  returns  position 
estimates  to  the  network.  In  a  distributed  implementation,  the  LaSLAT 
computations  are  performed  in-network. 

Consider  a  small  network  in  a  controlled  environment  such  as  a 
building.  The  experiments  in  chapter  6  are  all  representative  of  such  a 
scenario.  In  this  case,  the  central  computer  can  be  positioned  within 
a  single  radio  hop  of  all  or  nearly  all  sensors.  In  this  case,  the  best 
performance  is  obtained  by  transferring  all  the  measurements  to  the 
central  computer.  Centralization  in  this  case  reduces  the  radio  band¬ 
width,  memory,  and  processing  requirements  on  the  sensors,  decreasing 
the  hardware  cost  of  the  network.  Currently,  sensors  remain  somewhat 
expensive  and  radio  bandwidth  remains  a  scarce  commodity,  so  cen¬ 
tralization  is  very  desirable. 

As  the  network  grows  larger  distributed  computation  begins  to  look 
compelling.  When  multiple  hops  are  required  to  transmit  data  to  the 
central  computer,  the  energy  and  bandwidth  cost  of  centralization  in¬ 
creases,  particularly  for  nodes  near  the  central  computer.  In  this  case, 
distributing  LaSLAT  may  save  power  by  keeping  computation  and  com¬ 
munication  local.  Furthermore,  in  some  scenarios  (such  as  battlefield 
applications),  it  may  be  impossible  to  provision  the  network  with  a 
central  computer,  in  which  case  distributed  computation  is  necessary. 

In  many  cases,  however,  centralized  computation  is  straightforward 
and  desirable  even  in  large  networks.  For  instance,  in  large  indoor 
networks,  it  is  possible  to  use  pre-existing  high-speed  wireless  or  wired 
networks  to  transmit  data  rapidly  to  a  central  computer,  which  need 
not  be  co-located  with  the  network.  In  these  environments,  the  network 
forms  a  hierarchy  in  which  sensors  transmit  directly  to  a  base  station, 
which  in  turn  transmits  to  the  central  computer  over  a  high-speed  link. 
This  architecture  has  the  advantage  that  the  individual  sensors  need 
relatively  few  capabilities.  The  tracking  sensor  required  for  LaSLAT, 
for  example,  could  be  little  more  than  a  radio,  an  ultrasound  receiver, 
and  a  tone  detector.  Such  a  limited  sensor  is  likely  to  be  cheaper 
and  require  less  power  than  a  sensor  with  enough  processing  power, 
memory,  and  radio  bandwidth  to  perform  distributed  computations. 
The  savings  may  be  used  to  provision  a  smaller  number  of  computation 
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nodes  or  base  stations  with  a  high  speed  link  or  fast  processor  and  a 
more  substantial  power  supply.  Furthermore,  the  hierarchical  network 
is  better  suited  to  delivering  tracking  data  from  the  network  to  end 
users,  since  end  users  are  more  likely  to  be  connected  to  a  high  speed 
network  than  to  the  sensor  network’s  radio  system. 

LaSLAT  has  been  designed  to  perform  well  with  either  a  central¬ 
ized  or  a  distributed  implementation.  The  results  in  this  thesis  were  all 
computed  using  the  centralized  variant.  The  distributed  variant  per¬ 
forms  exactly  the  same  computations,  and  can  therefore  be  expected 
to  achieve  the  same  results. 
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Chapter  6 

Results 


Our  experiments  use  the  Cricket  ranging  system  [1]  (see  figure  6.1). 
Sensor  Crickets  are  placed  in  an  area,  and  one  Cricket  is  attached  to 
a  mobile.  The  mobile  Cricket  periodically  emits  an  event  (a  radio  and 
ultra-sound  pulse)  at  a  rate  between  one  and  three  per  second.  At  each 
sensor,  the  difference  in  arrival  time  of  these  two  signals  is  propor¬ 
tional  to  the  distance  between  the  sensor  and  mobile.  The  crickets  can 
therefore  estimate  ranges  from  these  arrival  times.  No  range  measure¬ 
ments  between  the  sensor  Crickets  are  collected.  The  measurements 
are  transmitted  to  a  desktop  machine,  which  processes  them  in  batches 
using  LaSLAT,  which  we  implemented  in  Java.  The  ultra-sound  sensor 
on  a  Cricket  occupies  a  1  cm  by  2  cm  area  on  the  circuit  board,  so  it 


Figure  6.1:  The  Cricket  sensor  node  used  in  our  experiments. 
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Figure  6.2:  Small  network  setup.  Six  sensors  (squares)  are  arranged 
around  a  rectangular  enclosure.  A  camera  captured  the  ground  truth 
trajectory  of  the  ROOMBA.  The  ROOMBA  followed  the  trajectory 
depicted. 


difficult  to  estimate  the  ground  truth  location  of  a  Cricket  beyond  that 
accuracy. 

In  all  of  the  experiments,  LaSLAT  recovered  sensor  locations  in 
a  relative  coordinate  system  that  can  be  aligned  to  a  global  coordi¬ 
nate  system  using  a  single  rigid  transformation  (a  rotation  and  transla¬ 
tion).  In  order  to  compare  LaSLAT’s  results  to  measured  ground  truth, 
we  computed  the  necessary  rigid  transformation  using  the  method  de¬ 
scribed  in  [15]. 

6.1  ROOMBA  Experiment 

Our  first  experiment  used  the  same  setup  as  [22].  Six  sensor  crickets 
were  placed  around  a  rectangular  enclosure  2.1  meters  by  1.6  meters. 
A  ROOMBA  robotic  vacuum  cleaner  with  the  mobile  Cricket  attached 
was  allowed  to  move  freely  within  the  enclosure,  generating  about  250 
events.  See  Figure  6.2. 

Most  events  were  measured  by  all  6  sensors.  An  initial  localization 
guess  was  obtained  from  radio  connectivity  information  using  the  ini- 
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Figure  6.3:  Recovered  trajectory  and  sensor  positions.  Circles  are 
guesses  of  initial  sensor  locations  obtained  from  radio  connectivity  (ap¬ 
pendix  A).  LaSLAT  processed  measurements  in  batches  of  30  events, 
and  recovered  sensor  locations  depicted  by  crosses.  The  trajectory  is 
also  correctly  recovered.  LaSLAT  improves  considerably  on  the  initial 
localization  guess  obtained  from  connectivity.  After  a  global  rotation 
and  translation,  the  average  localization  error  for  the  sensors  was  1.8 
cm,  which  is  within  the  error  tolerance  of  the  ground  truth. 
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tialization  routine  of  [29].  The  resulting  average  localization  error  in 
this  initial  guess  was  66  cm.  This  initial  guess  was  used  as  a  prior  and  an 
initial  iterate  for  LaSLAT.  LaSLAT  incorporated  range  measurements 
in  batches  of  30.  Each  mode  finding  operation  required  an  average  of 
only  2.8  Newton-Raphson  steps.  Figure  6.3  shows  the  estimated  sensor 
localizations  and  trajectory.  Since  the  output  had  an  arbitrary  rota¬ 
tion  and  translation,  it  was  rigidly  aligned  to  fit  the  rotation  and  origin 
of  the  ground  truth  using  the  algorithm  described  in  [15].  The  final 
localization  error  was  1.8  cm,  averaged  over  the  sensor  nodes.  This  is 
within  the  error  tolerance  of  the  ground  truth. 


6.2  Two  Dimensional  Experiments 

In  the  remainder  of  our  two  dimensional  experiments,  sensors  were 
placed  facing  upwards  on  the  floor  of  the  coverage  area.  The  mobile 
cricket  was  manually  moved  through  the  network.  It  was  suspended 
at  a  constant  height  of  about  190  centimeters,  and  oriented  facing  the 
floor.  Due  to  the  conical  propagation  of  ultrasound  from  the  mobile, 
this  approximated  a  radial  spread  of  sound  in  two  dimensions.  Range 
measurements  gathered  in  these  experiments  were  adjusted  in  a  pre¬ 
processing  phase  to  remove  the  effect  of  relative  height. 

6.2.1  27  node  experiment 

Our  second  experiment  involved  a  larger  network  with  27  Cricket  sen¬ 
sor  nodes  deployed  in  a  7  m  by  7  m  room.  Whereas  in  the  previous 
experiment  the  nodes  were  on  the  perimeter  of  the  ROOMBA  mobile’s 
trajectory,  in  this  experiment,  we  manually  pushed  a  mobile  through 
the  network,  generating  about  1500  events.  Figure  6.4  shows  the  loca¬ 
tion  of  the  sensors  and  part  of  the  trajectory  of  the  mobile  projected 
on  a  top  view  picture  of  the  setup. 

Each  event  was  heard  by  about  10  sensors.  Figures  6.5(a)-(d)  show 
localization  and  tracking  output  as  event  batches  are  processed,  along 
with  the  ground  truth  and  estimated  mobile  trajectories  for  that  batch. 
Error  ellipses  show  unit  standard  deviation  contours  for  each  sensor 
node.  Nodes  have  high  uncertainty  at  early  stages,  but  when  the  mo¬ 
bile  passes  near  a  node,  its  error  ellipse  shrinks  appropriately.  In  this 
experiment,  a  measurement  bias  of  about  23  cm  was  computed  for  each 
sensor  node.  Figure  6.6  shows  the  final  localization  of  the  nodes,  repro¬ 
jected  on  the  picture  of  the  setup.  This  experiment  used  batches  of  10 
measurements  and  produced  a  final  localization  error  of  1.9  cm.  Since 
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Figure  6.4:  Sensor  locations  and  mobile  trajectory  for  a  medium  size 
network.  Circles  outline  each  of  the  27  sensor  nodes.  Markers  on  the 
trajectory  depict  the  location  of  events.  250  of  the  1500  events  are 
shown,  with  consecutive  events  connected  by  a  line.  The  mobile  was 
offset  from  the  ground  plane  and  could  pass  over  nodes.  To  generate 
this  figure,  a  homography  that  accounts  for  the  camera  transformation 
was  used  to  project  real-world  coordinates  to  image  coordinates. 
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Figure  6.5:  The  output  of  LaSLAT  after  incorporating  (a)  50,  (b)  120, 
(c)  160,  and  (d)  1510  events.  The  batch  size  was  10.  Recovered  mo¬ 
bile  trajectory  (crosses)  and  ground  truth  mobile  trajectory  (solid  line) 
for  the  latest  batch  are  connected  by  a  line  to  show  correspondences. 
Estimated  mobile  locations  (dark  rectangles)  and  the  ground  truth  mo¬ 
bile  locations  (light  rectangles)  are  also  connected  with  a  line  to  show 
correspondence.  Error  ellipses  shrink  as  more  data  becomes  available. 
Between  events  120  and  160  (sub-figures  (b)  and  (c)),  the  mobile  swept 
around  the  bottom  of  the  network,  and  the  error  ellipses  and  localiza¬ 
tion  error  diminished  for  those  sensors.  Tracking  improved  as  sensors 
became  better  localized. 
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Figure  6.6:  Final  LaSLAT  localization  result,  with  batch  size  of  10. 
Crosses  show  estimated  sensor  locations.  These  are  correctly  estimated 
to  fall  on  the  corresponding  sensor.  Average  localization  is  1.9  cm. 


the  ground  truth  is  only  accurate  to  a  few  centimeters,  localization 
performance  is  best  examined  visually  via  Figure  6.6. 

We  compared  LaSLAT  using  varying  batch  sizes  to  the  Extended 
Kalman  Filter  (EKF),  which  is  identical  to  LaSLAT  limited  to  one 
Newton-Raphson  iteration.  Figure  6.7  shows  average  localization  er¬ 
rors  as  events  were  processed.  The  EKF  performs  best  with  no  batch¬ 
ing  (batch  size  =  1).  LaSLAT  converges  faster  and  also  exhibits  lower 
steady  state  localization  error.  As  batch  sizes  are  increased,  so  does 
the  rate  of  convergence  of  LaSLAT.  Batching  also  improves  the  final 
localization  error.  LaSLAT,  with  batch  sizes  of  1,  10  and  40,  produced 
final  localization  errors  of  3  cm,  1.9  cm,  and  1.6  cm  respectively.  On  av¬ 
erage,  LaSLAT  took  3  Newton-Raphson  iterations  to  incorporate  each 
batch.  The  EKF’s  final  localization  error  was  7.5  cm,  which  is  outside 
the  error  tolerance  for  the  ground  truth. 
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Figure  6.7:  Localization  error  as  a  function  of  the  number  of  events  ob¬ 
served  for  EKF  and  various  batch  sizes  for  LaSLAT.  LaSLAT  converges  more 
quickly  and  attains  a  lower  steady  state  error  than  the  EKF.  Furthermore, 
larger  batch  sizes  improve  the  convergence  rate  and  the  steady  state  error  of 
LaSLAT. 


44 


Figure  6.8:  A  sparser  2D  sensor  network  with  49  nodes  in  a  10m  by 
17m  environment. 
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Figure  6.9:  LaSLAT  localization  result  on  a  sparser  sensor  network  with 
49  nodes  in  a  10m  by  17m  environment.  Crosses  indicate  the  recovered 
sensor  locations,  projected  onto  the  image.  The  average  localization 
error  was  7.5  cm. 
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6.2.2  49  node  experiment 

Figure  6.9  shows  localization  results  on  a  larger  network  (49  sensors) 
deployed  over  a  larger  area  (10  m  by  17  nr).  For  comparison,  the 
experiment  setup  is  shown  in  figure  6.8.  With  about  0.3  sensors  per 
square  meter,  this  network  is  about  half  as  dense  as  the  one  shown 
in  Figure  6.4,  which  had  about  0.5  sensors  per  square  meter.  As  a 
result,  on  average  only  5  sensors  heard  each  event,  and  the  localization 
error  was  about  7.5  cm.  The  algorithm  also  determined  a  measurement 
biases  of  about  20  cm  for  all  nodes.  For  all  batch  sizes,  the  EKF 
produced  an  average  localization  error  of  about  80  cm,  showing  that 
the  improvement  due  to  Laplace’s  method  can  be  very  significant. 


6.3  Three  Dimensional  Experiments 

It  is  possible  to  use  LaSLAT  to  localize  and  track  sensors  and  mobiles 
in  a  three  dimensional  environment.  We  performed  two  experiments 
designed  to  test  LaSLAT  in  reasonable  three  dimensional  deployments. 

For  these  deployments,  crickets  were  placed  on  the  floor  as  in  the 
two  dimensional  experiments.  In  addition,  crickets  were  attached  to 
walls  using  Velcro.  A  few  sensors  were  also  placed  atop  furniture  in  the 
area.  We  also  modified  the  mobile  cricket  by  attaching  two  additional 
ultrasound  transducers.  These  transducers  more  closely  simulate  an 
omnidirectional  acoustic  pulse  than  the  conic  emanation  of  the  stan¬ 
dard  cricket  transducer.  See  figure  6.10  for  a  picture  of  this  modified 
cricket. 

The  three  dimensional  environments  present  two  obstacles  that  are 
not  present  in  the  two  dimensional  experiments.  First,  the  furniture  in 
the  area  drastically  increases  the  effect  of  echoes  on  range  measurement 
quality.  These  echoes  are  suppressed  using  the  measurement  outlier  re¬ 
jection  scheme  presented  in  section  4.1.  Second,  the  initialization  algo¬ 
rithm  we  used  in  two  dimensions  [29]  does  not  provide  for  localization 
in  3D.  However,  we  developed  a  simple  extension  (see  appendix  A)  to 
allow  us  to  initialize  LaSLAT  satisfactorily. 

6.3.1  40  node  experiment 

In  our  first  three  dimensional  experiment,  we  placed  40  crickets  on  the 
floor  and  walls  of  a  4  x  6  meter  room,  which  contained  all  of  its  normal 
furniture:  tables,  chairs,  printers,  and  a  refrigerator.  This  mobile  was 
carried  by  hand  through  the  room  and  moved  completely  arbitrarily, 
including  changes  in  speed,  loops,  and  twists.  Each  event  was  observed 
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Figure  6.11:  Smaller  environment  used  for  3D  localization  experiments. 
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Figure  6.12:  LaSLAT  results  plotted  on  a  picture  of  the  network.  Plus 
signs  indicate  the  estimated  3D  positions  of  the  sensors.  A  small  por¬ 
tion  of  the  mobile  trajectory  (about  80  events)  is  plotted  as  asterisks 
connected  by  a  dotted  line.  LaSLAT  localized  sensors  to  within  7  cm. 


by  on  average  17  sensors.  The  network  is  shown  in  figure  6.11.  Lo¬ 
calization  results  are  plotted  in  figure  6.12.  LaSLAT  localized  sensors 
to  within  7  cm  while  successfully  tracking  the  path  of  the  mobile  in 
3D.  Much  of  this  error  is  accounted  for  by  the  difficulty  of  measuring 
ground  truth  in  this  environment. 

The  best  3D  results  were  obtained  using  a  relatively  large  batch 
size  of  250  events.  Smaller  batch  sizes  lose  too  much  information  in 
the  Gaussian  approximation  used  to  preserve  state  between  batches, 
causing  LaSLAT  to  converge  slowly.  In  3D,  the  optimization  (3.6)  has  a 
more  complicated  structure  due  to  the  additional  unknown  parameters, 
which  larger  batches  help  clarify.  Thus,  a  large  initial  batch  helps 
LaSLAT  make  rapid  progress  towards  an  accurate  estimate.  After  the 
first  few  batches,  the  batch  size  may  be  decreased  to  speed  computation, 
since  the  later  batches  need  only  refine  the  already  fairly  high  quality 
estimate. 

The  3D  environment  shown  in  figure  6.11  is  considerably  less  san¬ 
guine  to  ultrasound  time-of-flight  ranging  than  the  two  dimensional  en¬ 
vironments.  Most  events  caused  several  very  erroneous  measurements. 
Figure  6.13  shows  a  partial  screen  capture  from  our  Java-based  LaSLAT 
implementation.  The  dark  line  indicates  a  portion  of  the  trajectory  re¬ 
covered  by  LaSLAT.  The  light  line  is  recovered  using  multilateration  on 
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Figure  6.13:  Sample  LaSLAT  3D  tracking  results.  The  black  line  shows 
the  mobile  trajectory  recovered  by  LaSLAT  on  the  40  node  3D  topology. 
The  gray  line  shows  the  mobile  trajectory  recovered  using  multilater- 
ation  on  raw  measurements  and  ground  truth  sensor  positions.  Stars 
represent  the  position  of  the  mobile  at  the  time  of  an  event.  Approx¬ 
imately  three  events  occurred  every  second.  Note  that  LaSLAT  very 
noticeably  outperforms  the  naive  tracking  algorithm,  even  though  it 
does  not  know  the  sensor  positions  a  priori. 
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the  ground  truth  sensor  positions  and  the  raw  measurement  data,  and 
therefore  serves  as  a  baseline  tracking  algorithm.  LaSLAT’s  improved 
tracking  performance  is  due  to  the  combination  of  sensor  calibration 
(section  3.2),  measurement  outlier  rejection  (section  4.1),  and  smooth 
dynamics  (section  4.2).  Some  of  these  techniques  can  be  and  have  been 
adapted  to  pure  tracking  applications.  However,  LaSLAT  is  able  to 
perform  the  same  high  fidelity  tracking  while  simultaneously  localizing 
and  calibrating  the  network. 

We  ran  LaSLAT  both  with  and  without  smooth  dynamics  and  out¬ 
lier  rejection  to  observe  the  techniques’  impacts  on  localization  accu¬ 
racy.  The  results  for  the  40-node  3D  experiment  are  graphed  in  figure 
6.14.  Note  that  LaSLAT  with  smooth  dynamics  and  outlier  rejection 
performs  best,  with  a  final  accuracy  of  less  than  7  centimeters. 

Smooth  dynamics  makes  relatively  little  difference  to  sensor  local¬ 
ization,  as  its  final  average  position  error  is  8  centimeters.  This  is  ac¬ 
counted  for  by  the  fact  that  smooth  dynamics  primarily  affects  events 
that  deviate  from  their  neighbors.  Since  most  events  tend  to  be  ac¬ 
curately  placed  without  requiring  smooth  mobile  dynamics,  correcting 
the  few  deviants  events  has  only  a  small  effect  on  the  sensor  position 
estimates. 

However,  measurement  outlier  rejection  makes  a  very  substantial 
difference  on  this  data  set.  Without  outlier  rejection,  the  final  localiza¬ 
tion  estimate  has  a  mean  error  of  31  centimeters.  This  suggests  that 
the  outlier  rejection  technique  presented  in  section  4.1  is  effective  at 
detecting  and  de-emphasizing  errant  measurements. 

For  reference,  we  also  plot  the  performance  of  the  Extended  Kalman 
Filter  (EKF)  on  this  40-node  experiment.  The  final  EKF  estimate  has 
on  average  over  half  a  meter  of  localization  error. 

6.3.2  55  node  experiment 

In  our  second  3D  experiment,  we  used  55  sensors  to  cover  the  floor 
and  walls  of  a  7  x  10  meter  room.  Since  the  sensors  were  packed  less 
densely,  on  average  only  11  sensors  witnessed  each  event.  A  portion  of 
the  network  is  shown  in  figure  6.15,  and  the  results  are  plotted  in  figure 
6.16.  LaSLAT  localized  the  sensors  to  within  7  cm,  consistent  with  the 
previous  experiment. 

This  environment  was  better  suited  to  the  Cricket  ultrasound  rang¬ 
ing  system  since  it  contained  few  obstacles.  Consequently,  sensors  re¬ 
ported  fewer  bad  measurements.  Nevertheless,  enough  poor  measure¬ 
ments  were  taken  that  outlier  rejection  improved  performance  notice¬ 
ably. 
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Figure  6.14:  Performance  impact  of  dynamics  and  outlier  rejection  on 
a  3D  dataset.  Note  that  all  four  variants  of  LaSLAT  outperform  the 
Extended  Kalman  Filter  (EKF)  implementation. 


Figure  6.15:  Part  of  the  55  node  network  for  3D  localization  experi¬ 
ments.  Sensors  are  circled  in  blue. 
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Figure  6.16:  Results  for  the  55  node  3D  experiment.  Estimated  posi¬ 
tions  are  plotted  with  crosses,  and  a  portion  of  the  recovered  mobile 
trajectory  is  plotted  using  asterisks  and  a  dotted  line.  LaSLAT  local¬ 
ized  this  network  with  a  mean  error  of  7  centimeters. 
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Chapter  7 

Future  Work 


In  addition  to  the  features  developed  in  this  thesis,  LaSLAT  raises 
several  interesting  possibilities  for  future  work. 


7.1  Hallway  Alignment 

As  we  demonstrated  in  chapter  6,  LaSLAT  perforins  very  well  in  net¬ 
works  that  have  a  reasonably  convex  shape.  However,  its  performance 
diminishes  when  the  network  is  divided  into  convex  areas  separated 
by  narrow  connections:  for  instance,  two  rooms  separated  by  a  nar¬ 
row  hallway.  In  these  topologies,  LaSLAT  accurately  localizes  sensors 
in  the  rooms  with  respect  to  other  sensors  in  the  same  room.  How¬ 
ever,  it  cannot  always  orient  the  rooms  as  a  whole  with  respect  to  each 
other.  It  appears  that  the  measurements  gathered  from  events  in  and 
around  the  hallway  area  typically  fail  to  unambiguously  resolve  the 
alignment.  Augmenting  LaSLAT  using  a  coordinate  system  stitching 
approach  similar  to  [7,21,22,24,31]  could  help  resolve  this  ambiguity. 

7.2  Stationary  Targets 

In  our  experiments,  the  mobile  moved  continuously  for  the  majority 
of  the  experiment.  However,  real  mobiles  can  be  expected  to  remain 
stationary,  potentially  for  lengthy  periods  of  time.  After  LaSLAT  has 
been  running  for  some  time,  this  poses  no  problem.  However,  during 
the  first  few  batches  of  LaSLAT  difficulties  can  arise.  A  stationary 
mobile  causes  a  large  number  of  highly  similar  events.  This  results  in 
a  batch  of  events  that  resolves  few  ambiguities  in  the  sensor  positions. 
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This  leads  in  turn  to  a  multi-modal  posterior  whose  complexity  is  lost 
in  LaSLAT’s  Gaussian  approximation.  As  a  result,  LaSLAT’s  estimate 
may  converge  artificially  to  a  very  poor  estimate.  The  corresponding 
prediction  has  low  covariance  since  it  is  the  product  of  many  events, 
but  is  nevertheless  quite  inaccurate. 

There  are  a  number  of  techniques  available  to  resolve  this  prob¬ 
lem.  First,  it  may  be  possible  to  use  a  heuristic  over  the  collected 
measurements  or  the  LaSLAT  results  to  detect  uninformative  batches 
of  events.  Bad  batches  can  then  be  ignored.  In  effect,  LaSLAT  would 
wait  to  begin  localization  until  good  quality  data  becomes  available. 
Alternatively,  we  could  use  a  multi-modal  posterior  approximation  in¬ 
stead  of  the  Gaussian  approach  detailed  in  this  thesis. 

7.3  Multi-modal  Posterior  Approximation 

The  Gaussian  approximation  described  in  section  3.3  has  an  important 
disadvantage:  it  cannot  accurately  represent  a  multi-modal  posterior 
distribution.  Consequently,  LaSLAT  requires  that  the  initial  measure¬ 
ment  batches  identify  a  single  high  probability  mode.  This  necessity 
can  be  diminished  by  using  a  different  approximation  such  as  a  mix¬ 
ture  of  Gaussians  or  a  particle  filter.  These  approximations  can  model 
multi-modal  distributions;  however,  they  increase  the  processing  load 
for  a  LaSLAT  batch. 

We  were  not  able  to  experiment  with  these  types  of  representation 
for  this  thesis.  However,  they  may  be  necessary  for  accurate  localization 
of  some  networks,  especially  when  areas  have  few  measurements  or  the 
mobile  is  somewhat  stationary.  They  may  also  permit  reduced  batch 
sizes  when  performing  three  dimensional  SLAT. 


7.4  Low  Quality  Hardware 

In  this  thesis,  we  used  the  Cricket  platform  [1]  for  all  of  our  experi¬ 
ments.  The  Cricket  ranging  algorithm  is  highly  accurate  in  the  absence 
of  echoes,  achieving  1-2  centimeter  accuracy  on  ranges  of  up  to  ten  me¬ 
ters.  However,  we  suspect  that  LaSLAT  will  remain  effective  when  less 
accurate  hardware  is  used. 

One  significant  cause  of  error  in  some  acoustic  time  of  flight  rang¬ 
ing  systems  can  be  manufacturing  differences  between  sensors.  In  our 
experiments  with  the  Cricket  system,  we  encountered  no  need  to  model 
such  inter-sensor  differences  explicitly.  In  fact,  as  we  indicated  in  sec¬ 
tion  3.5,  we  found  it  helpful  to  encourage  sensor  measurement  biases 
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to  be  similar.  However,  [35]  shows  that  some  types  of  acoustic  ranging 
can  benefit  from  careful  modeling  of  individual  sensor  characteristics. 
LaSLAT  seems  well-suited  to  perform  this  type  of  analysis. 

We  have  shown  in  this  thesis  that  LaSLAT  performs  very  well  on 
Cricket  networks  that  are  deployed  with  reasonable  density  in  indoor 
environments;  however,  we  are  very  interested  in  determining  the  min¬ 
imal  hardware  profile  required  for  effective  use  of  LaSLAT.  This  means 
performing  additional  experiments  on  networks  with  reduced  sensor 
density,  event  frequency,  and  measurement  accuracy.  The  low  density 
experiments  presented  in  chapter  6  represent  initial  research  in  this 
area.  Nevertheless,  much  remains  to  be  done,  especially  with  regard  to 
the  use  of  cheaper  and  lower  quality  sensor  hardware. 


7.5  Sensor  Dynamics 

LaSLAT  is  well-suited  to  localizing  moving  sensors.  In  section  3.4,  we 
stated  the  assumption  that  sensor  positions  are  stationary.  However, 
LaSLAT  could  easily  be  adapted  to  use  the  more  relaxed  constraint 
that  sensors  may  move  over  time.  This  sensor  mobility  constraint  might 
resemble  the  smooth  dynamics  extension  discussed  in  section  4.2.  This 
would  allow  LaSLAT  based  tracking  networks  even  greater  flexibility: 
a  user  could  move  sensors  after  deployment  to  optimize  the  network.  If 
sensors  were  able  to  move  autonomously,  sensors  could  automatically 
migrate  to  achieve  the  desired  quality  of  sensor  coverage. 

7.6  Acoustic-only  LaSLAT 

This  thesis  has  been  written  under  the  assumption  that  sensors  com¬ 
pute  ranges  using  Time  Difference  of  Arrival  (TDoA)  ranging;  thus, 
the  mobile  emits  a  tagged  radio  pulse  followed  by  an  acoustic  pulse. 
Sensors  use  the  difference  in  time  of  arrival  between  the  radio  message 
and  the  acoustic  pulse  to  estimate  distance. 

However,  the  measurement  model  (section  3.2)  can  be  adapted  to 
use  an  acoustic  pulse  alone.  In  this  case,  the  sensors  must  be  time 
synchronized,  which  can  be  done  using  an  off-the-shelf  time  synchro¬ 
nization  system  such  as  [18].  The  sensors  measure  the  arrival  times 
of  acoustic  pulses,  which  they  report  to  LaSLAT.  These  arrival  times 
can  be  seen  as  a  function  of  the  true  distance  between  the  event  and 
the  sensor  and  the  true  time  of  the  event.  The  event  time  becomes  an 
additional  parameter  to  be  estimated  by  the  Bayesian  filter.  Such  a 
measurement  model  would  allow  SLAT  to  localize  and  track  using  only 
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acoustic  time  of  arrival.  More  generally,  it  may  be  possible  to  localize 
using  ambient  environmental  noise,  as  long  as  recognizable  noises  are 
witnessed  by  enough  sensors. 


7.7  Distributed  Implementation 

The  results  provided  in  chapter  6  were  generated  using  a  centralized 
implementation.  Measurements  were  taken  in  a  real  network,  then  ex- 
filtrated  to  a  laptop  for  processing.  However,  as  we  show  in  chapter 
5,  LaSLAT  can  be  feasibly  distributed  to  perform  the  necessary  com¬ 
putations  in-network.  As  we  argue  in  that  section,  it  is  not  always 
desirable  to  distribute  LaSLAT  processing;  however,  a  distributed  im¬ 
plementation  will  allow  LaSLAT  to  be  of  assistance  in  networks  where 
centralized  computation  is  infeasible. 
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Chapter  8 

Conclusion 


In  this  thesis  we  presented  LaSLAT,  a  sensor  network  algorithm  that  si¬ 
multaneously  localizes  sensors,  calibrates  sensing  hardware,  and  tracks 
unconstrained  moving  targets  using  only  range  measurements  between 
the  sensors  and  the  mobile.  On  both  two-  and  three-dimensional  net¬ 
works  using  the  Cricket  ranging  system,  LaSLAT  was  able  to  localize 
sensors  to  within  several  centimeters  of  their  ground  truth  positions 
while  recovering  a  range  measurement  bias  for  each  sensor  and  the 
complete  trajectory  of  the  target. 

LaSLAT  is  based  on  a  Bayesian  filter,  which  updates  a  probabil¬ 
ity  distribution  over  the  quantities  of  interest  as  measurements  arrive. 
The  algorithm  is  distributable,  and  requires  only  a  constant  amount 
of  space  with  respect  to  the  number  of  measurements  incorporated. 
LaSLAT  is  easy  to  adapt  to  new  types  of  hardware  and  new  physical 
environments  due  to  its  use  of  intuitive  probability  distributions:  one 
adaptation  demonstrated  in  this  thesis  uses  a  mixture  measurement 
model  to  detect  and  compensate  for  bad  acoustic  range  measurements 
due  to  echoes. 

Localization  and  calibration  using  an  unconstrained  mobile  offers 
several  advantages  over  localization  and  calibration  using  inter-sensor 
ranges.  Since  the  mobile  may  move  freely,  the  measurements  taken 
by  the  sensors  have  greater  diversity  and  therefore  diminished  likeli¬ 
hood  of  systematic  errors.  In  particular,  the  mobile  can  help  overcome 
interference  from  fixed  obstacles  such  as  furniture. 

Use  of  a  mobile  also  means  that  the  sensors  are  no  longer  con¬ 
strained  to  have  line  of  sight  to  each  other.  This  enables  a  broader 
range  of  network  deployments:  in  particular,  it  allows  sensors  to  be 
placed  to  optimize  tracking  coverage  rather  than  to  accommodate  the 
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peculiarities  of  a  localization  algorithm. 

LaSLAT  in  its  current  form  could  be  used  to  facilitate  deployment  of 
indoor  navigation  systems,  either  for  human  use  or  as  a  guidance  system 
for  mobile  robots.  Using  LaSLAT,  the  sensors  for  such  a  system  could 
be  deployed  in  an  ad-hoc  fashion  on  walls,  floors,  or  ceilings.  They 
could  then  use  the  movements  of  persons  to  be  tracked  (i.e.  targets)  to 
localize  and  calibrate  in-place,  forming  a  high-precision  location  service 
with  little  human  effort. 

LaSLAT  also  suggests  an  interesting  automated  sensor  deployment 
system  using  a  mobile  robot.  Such  a  robot  would  survey  its  environ¬ 
ment,  periodically  dropping  sensors.  The  sensors  would  provide  high 
quality  landmarks  for  the  robot’s  mapping  system,  while  the  robot  pro¬ 
vided  a  moving  target  for  LaSLAT-based  localization  of  the  sensors. 
The  sensors  could  even  provide  a  long  term  record  of  the  mapping  op¬ 
eration,  which  they  could  propagate  to  other  robots  not  tasked  with 
surveying. 

While  many  impediments  remain,  LaSLAT  is  an  important  step  to¬ 
ward  cheap,  easily  deployable,  and  accurate  sensor  networks  for  target 
tracking. 
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Appendix  A 

Initialization 


Initialization  plays  an  important  role  in  LaSLAT.  As  several  authors 
[22, 29]  have  observed,  inferring  locations  using  distance  measurements 
can  suffer  from  local  minima.  The  Newton-Raphson  optimization  de¬ 
scribed  in  appendix  B  does  not  always  find  the  globally  minimal  posi¬ 
tioning.  Rather,  it  finds  a  minimum  near  its  initial  value,  which  may 
or  may  not  be  the  global  minimum. 

This  problem  is  typically  solved  in  one  of  two  ways.  The  first 
method  involves  running  the  optimization  several  times  using  different 
random  choices  of  initial  value.  The  best  result  from  these  optimiza¬ 
tions  is  taken  as  the  global  minimum. 

In  the  second  method,  which  we  use  in  LaSLAT,  the  initial  iterate 
is  chosen  intelligently  to  be  close  to  the  true  solution.  As  a  result,  the 
Newton-Raphson  optimization  converges  to  the  global  minimum  on  the 
first  try.  After  the  first  LaSLAT  iteration,  the  estimate  from  the  pre¬ 
vious  iteration  is  generally  a  good  initial  iterate  for  Newton-Raphson. 
However,  we  depend  on  a  separate  initialization  algorithm  to  provide 
an  initial  iterate  for  the  first  LaSLAT  iteration.  In  particular,  we  use  a 
distributed  initialization  algorithm  suggested  by  Priyantha  et  al.  [29]. 
We  have  also  experimented  with  techniques  such  as  gradient  multi- 
lateration  [23]  and  multidimensional  scaling  [17].  In  our  experience, 
the  algorithm  from  [29]  offered  the  best  performance.  This  algorithm 
provides  a  coarse  position  estimation  for  the  sensors  in  the  network 
by  defining  a  polar  coordinate  system  based  on  hop  counts  between 
specially  chosen  sensors.  It  is  described  in  detail  below. 


60 


A.l  Initialization  Using  Radio  Connectiv¬ 
ity 

Let  hij  be  the  minimum  number  of  radio  hops  required  to  forward  a 
message  from  sensor  i  to  sensor  j.  Perform  the  following  steps: 

1.  Choose  node  no  to  be  the  sensor  with  minimal  ID.  Propagate  a 
gradient  from  no  to  compute  ho,k  for  all  k. 

2.  Choose  node  n\  to  be  the  sensor  k  with  maximal  break¬ 
ing  ties  by  ID.  m  is  one  “corner”  of  the  network.  Propagate  a 
gradient  from  n\. 

3.  Choose  node  n2  to  be  the  sensor  k  with  maximal  hi  k,  breaking 
ties  by  ID.  n,2  is  the  corner  opposite  n\.  Propagate  a  gradient 
from  7i2- 

4.  Choose  node  n 3  to  be  the  sensor  k  that  minimizes  | /i2,fc  —  /ii,fc|, 
breaking  ties  by  the  sum  h2,fc  +  (any  further  ties  are  resolved 
by  ID).  713  is  the  corner  between  n\  and  n2-  Propagate  a  gradient 
from  713. 

5.  Choose  node  714  to  be  the  sensor  k  that  minimizes  | h2,fc  —  /ii  fc|, 
breaking  ties  by  maximizing  hotk  (any  further  ties  are  resolved  by 
ID).  714  is  the  final  corner.  Propagate  a  gradient  from  714. 

6.  Choose  node  715  to  be  the  sensor  k  that  minimizes  \h2ik  —  hitk\  + 
\hi,k  —  hz,k\,  breaking  ties  by  ID.  This  sensor  is  closest  to  the 
center  of  the  network.  Propagate  a  gradient  from  715. 

7.  For  each  sensor  k,  define  pk  =  h^,k  *  R ,  and  tan(0fc)  =  > 

where  R  is  the  maximum  range  of  the  sensor’s  radio.  Compute 
the  position  sjl  by  converting  the  polar  coordinate  (pk,0k)  to 
rectangular  coordinates. 

This  initialization  algorithm  tends  to  be  highly  inaccurate  (see  fig¬ 
ure  A.l).  However,  the  chosen  positions  tend  to  help  the  LaSLAT 
optimization  avoid  local  minima  that  decrease  accuracy.  Our  specific 
choice  of  initialization  is  somewhat  irrelevant  -  any  reasonably  accurate 
localization  algorithm  could  be  substituted. 
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Figure  A.l:  Initialization  results  for  the  27  node  2D  experimental  topol¬ 
ogy.  Initial  position  estimates  are  marked  with  circles.  Lines  are  drawn 
from  these  estimates  to  the  corresponding  ground  truth  sensor  loca¬ 
tions.  The  average  position  estimate  error  is  70  cm. 
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A. 2  Initialization  in  3D 

The  algorithm  described  in  the  previous  section  is  not  designed  to  lo¬ 
calize  sensor  configurations  that  have  height  differences.  In  general, 
most  two  dimensional  localization  algorithms  extend  poorly  to  practi¬ 
cal  three  dimensional  topologies.  Most  two  dimensional  topologies  re¬ 
semble  a  “cloud”  of  sensors  -  sensors  are  evenly  scattered  throughout  a 
generally  convex  area.  In  three  dimensions,  this  is  typically  infeasible, 
especially  indoors.  Instead,  practical  deployments  look  more  like  figure 
6.11,  where  sensors  are  placed  mostly  on  walls  and  either  the  floor  or 
the  ceiling.  As  a  result,  naive  three  dimensional  extensions  to  2D  lo¬ 
calization  algorithms  such  as  [11,17,22,29,31]  tend  to  perform  badly 
on  real  networks. 

We  have  developed  a  simple  method  of  extending  these  algorithms 
to  three  dimensions,  using  a  single  additional  bit  of  information  that 
indicates  whether  a  sensor  is  placed  on  a  floor  or  a  wall.  This  informa¬ 
tion  could  be  obtained  by  a  relatively  cheap  and  imprecise  orientation 
sensor  on  each  node  or  by  manual  configuration.  Also,  it  is  probable 
that  in  an  indoor  environment  the  form  factor  of  the  sensors  will  vary 
between  wall  sensors  and  floor  sensors,  making  manual  configuration 
straightforward . 

If  the  radio  range  is  sufficiently  small  (this  may  be  accomplished  by 
temporarily  reducing  the  power  of  the  sensors’  radios),  a  3D  topology 
like  figure  6.11  will  appear  to  be  two  dimensional  from  the  perspective 
of  radio  connectivity,  where  wall  sensors  appear  at  the  outer  extremities 
of  the  topology.  In  this  case,  a  two  dimensional  localization  algorithm 
such  as  the  one  described  in  the  previous  section  may  be  applied.  As  a 
post-processing  step,  the  regions  of  the  two  dimensional  plane  contain¬ 
ing  wall  sensors  may  be  “folded”  upwards  into  the  height  dimension  to 
approximate  the  structure  of  the  3D  environment. 

In  our  experience  so  far,  this  trivial  modification  allows  algorithms 
such  as  [29]  to  be  successfully  used  in  3D  LaSLAT. 
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Appendix  B 

Newton- Raphson 


In  this  appendix,  we  show  how  to  optimize  equation  (3.6)  using  the 
Newton-Raphson  iterative  optimization  algorithm.  We  also  derive  the 
curvature  of  the  negative  log  posterior  required  for  Laplace’s  method 
(section  3.3),  and  show  that  mode  finding  preserves  local  connectivity. 


B.l  Finding  a  Mode 


Equation  (3.6)  can  be  expressed  in  the  form  of  non-linear  least  squares. 
Let  fij(x)  =  ||sf — e*-||+sf,  and  define  f(x)  as  a  column  vector  consisting 
of  all  fij  (x) .  Using  Ux  =  Cov-1  [x4|yoW],  and  /zx  =  E  [x*  |y°<rf] ,  we 
can  write  equation  (3.6)  as: 

argmin— ^ ||f(x)  -  y4 ||2  +  (x  -  /rx)Tnx(x  -  /xx)  (B.l) 
x  az 


Each  iteration  of  Newton-Raphson  maps  an  iterate  x ^  to  the  next 
iterate  x(4+b  by  approximating  (B.l)  with  a  linearization  about  xy>, 
then  optimizing  over  x: 


X(‘+D 


.  1 
arg  mm  — 

x  az 


Vf(4)x-  b 


+  (x  -  /rx)TUx(x  -  /xx), 


(B.2) 


where  the  matrix  Vfi^  is  the  derivative  of  f  with  respect  to  x  at  x^, 
and  the  column  vector  b  =  Vf^x^  —  f(x^4))  —  y4. 

Equation  (B.2)  is  a  linear  least  squares  problem  in  terms  of  x.  Its 
solution  can  be  found  by  setting  the  derivative  with  respect  to  x  to  zero. 
This  leaves  a  linear  problem  that  can  be  solved  by  matrix  inversion: 

[Ox  +  4  Vf{t)T  Vf(4)l  X  =  nxp  +  4  Vf  (t)Tb.  (B.3) 

c tz  crz 
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Furthermore,  differentiating  (B.2)  one  more  time  results  in  H  = 
flx  +  which  is  required  for  Laplace’s  method  (section 

3.3).  Since  (B.2)  is  an  approximation  to  the  negative  log  posterior 
(3.6),  H  serves  as  an  approximation  to  its  Hessian  at  x^b 


B.2  Locality  of  Mode  Finding 

Because  the  true  distance  /.i?  depends  only  on  sensor  i  and  event  loca¬ 
tion  j,  each  row  of  Vf^  is  made  up  of  zeros,  except  at  locations  cor¬ 
responding  to  the  zth  sensor  and  the  jth  event.  Thus  Vfd)TVf^)  has 
local  connectivity.  If  flx  =  Cov^1  [x*|yoW]  has  local  connectivity,  then 
the  updated  covariance  matrix  Cov-1  [x4|yoW]  =  flx  +  Vf^T Vf^/c2 
also  has  local  connectivity.  Therefore  as  claimed  in  chapter  5,  incorpo¬ 
rating  a  batch  of  measurements  preserves  local  connectivity. 
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